System and method for real time dynamic lighting simulation

ABSTRACT

Sustainable building lighting and energy modelling and control, and the associated computer graphics, including real-time dynamic lighting simulation, are concerned with: an optimized method for radiance modelling, including its application to predictive daylight harvesting; and the real-time simulation of physically-based electric lighting and daylighting for architectural, horticultural, and theatrical lighting systems visualization. In order to display and analyze in real time a photometrically accurate representation of an environment, thousands of lighting channels may have their intensity settings continually varied such that a user may interactively view the three-dimensional environment without the need for ongoing global illumination calculations. This can be accomplished utilizing texture maps as a multiplicity of canonical radiosity solutions, each representing a lighting channel for dynamic lighting simulation, and storing the solutions in the texture memory of a graphics processing unit.

This application claims benefit of U.S. provisional Ser. No. 62/369,912filed 2 Aug. 2016, which is incorporated by reference herein in itsentirety, and is a continuation in part of U.S. patent application Ser.No. 15/470,180 filed 27 Mar. 2017 which claims benefit of U.S.provisional Ser. No. 62/313,718 filed 26 Mar. 2016 and is a continuationin part of U.S. patent application Ser. No. 15/407,176 filed 16 Jan.2017 which claims benefit of U.S. provisional Ser. No. 62/279,764 filed17 Jan. 2016 and is a continuation in part of U.S. patent applicationSer. No. 14/792,590 filed 6 Jul. 2015 which claims benefit of U.S.provisional Ser. No. 62/172,641 filed 8 Jun. 2015 and is a continuationin part of U.S. patent application Ser. No. 13/446,577 filed 13 Apr.2012 (U.S. Pat. No. 9,078,299) which claims benefit of U.S. provisionalSer. No. 61/565,195 filed 30 Nov. 2011 and claims benefit of U.S.provisional Ser. No. 61/457,509 filed 14 Apr. 2011.

TECHNICAL FIELD

The subject matter of the present invention relates to the field ofsustainable building lighting and energy modelling and control andcomputer graphics, and more particularly, is concerned with: anoptimized method for radiance modelling, including its application topredictive daylight harvesting; and the real-time simulation ofphysically-based electric lighting and daylighting for architectural andtheatrical lighting systems visualization.

BACKGROUND

Electric lighting accounts for approximately 40 percent of all energyconsumed in modern buildings. Incorporating available daylight canreduce these annual energy costs by 40 to 50 percent using “daylightharvesting” techniques. The basic principle of daylight harvesting is tomonitor the amount of daylight entering an interior space and dim theelectric lighting as required to maintain a comfortable luminousenvironment for the occupants. Where required, motorized blinds andelectrochromic windows may also be employed to limit the amount ofdaylight entering the occupied spaces. Further energy savings can berealized through the use of occupancy sensors and personal lightingcontrols that operate in concert with the daylight harvesting system.These systems can be onerous to model and control using current methods,and rely on extensive computing resources to adequately model andcontrol the variables involved.

Three-dimensional virtual representations of architectural andtheatrical lighting systems were first developed in the 1980s, and todayare often used by architects and lighting designers to both visualizeand analyze their daylighting and electric lighting system designs.Architectural and theatrical lighting systems are not limited tointerior lighting, and include lighting systems on the exterior ofbuildings or structures.

A typical workflow for the creation of said virtual representationscomprises the steps of: 1) using a computer-aided drafting (CAD) programto create a three-dimensional representation of the geometry of abuilding or similar architectural environment; 2) assigning materialproperties such as color, transparency, and texture to the geometricsurfaces; 3) specifying the photometric properties of light sourceswithin the environment; 4) calculating the distribution of direct andindirect illumination due to the light sources within the environment;and 5) displaying the virtual representations on a computer monitor orother visual display medium.

Calculation of the distribution of direct and indirect illumination istypically performed using global illumination techniques such as thosedescribed in, for example, “The State of the Art in Interactive GlobalIllumination” (Ritschel et al., 2011, Computer Graphics Forum31(1):160-168). Critically for the photometric analysis of lightingdesigns (including both electric lighting and daylight), thesetechniques are based on the physical principles of light interactionwith surfaces and materials, rather than mere artistic representations.

Most global illumination techniques are based on ray tracing, includingthose techniques known as irradiance caching, photon mapping, andmultidimensional light cuts. The advantage of ray tracing techniques isthat they can accurately model reflections from specular andsemispecular (“glossy”) surfaces, refraction by transparent materialssuch as glass and water, interactions with participating media such assmoke and fog, and subsurface translucency of human skin, marble, and soforth. These features are critically important for the generation ofhigh quality “photorealistic” imagery, including both still images andmotion pictures.

One of the disadvantages of ray tracing techniques is that they areview-dependent; a specific position and view direction within thevirtual environment must be chosen in order to trace rays to and fromthe virtual camera. The various ray tracing techniques arecomputationally expensive, typically necessitating the distribution ofdirect and indirect illumination within the environment to be calculatedoffline using minutes to days of CPU time. (While it is possible toapproximate the appearance of direct and indirect illumination within anenvironment in real time for computer games, the resultant images areunsuitable for photometric analysis as required by architects andlighting designers.)

The view dependency of ray tracing techniques further means that onlythose surfaces visible in the rendered images are available forphotometric analysis. For an architectural building with multiple rooms,for example, one or more rendered images must be generated for eachroom.

Not all global illumination techniques are based on ray tracing; thereis also a class of finite element methods referred to as “radiosity.”These methods are more commonly used in thermal engineering to modelradiative flux transfer between surfaces, but they are equally capableof modeling reflections from diffuse and semispecular surfaces, andtransmission by transparent and translucent surfaces, as described in,for example, “Radiosity: A Programmer's Perspective” (Ashdown 1994. NewYork, N.Y.: John Wiley & Sons). While less capable of generating highquality photorealistic imagery, they are ideally suited for photometricanalysis. In addition to the calculated illuminance distributions beingphysically accurate, they are view-independent. Once the calculationshave been completed, the entire three-dimensional environment can beinteractively viewed and photometrically analyzed in real time withoutthe need for further calculations.

A particular disadvantage of both ray tracing techniques and radiositymethods, however, is that they are incapable of generatingrepresentations of “dynamic” lighting systems in real time. Examples ofsuch lighting systems include theatrical and entertainment lightingwhere the intensity and color of the light sources is continuouslyvaried, and architectural lighting systems where the light sources maybe dimmable or have varying color temperatures (e.g., “warm white” to“cool white”). Another example is daylighting, where it may be desirableto display and analyze the distribution of daylight in architecturalspaces throughout the year.

Theatrical and entertainment lighting systems in particular may havetens to hundreds, and even thousands, of lighting “channels” whereineach channel dims one or a plurality of electric light sources(“luminaires”). Three such channels may be assigned to color-changingluminaires with red, green, and blue light sources. Architecturallighting systems are similar, albeit with only a few lighting channelsto dim the luminaires and possibly change their color temperatures.

If possible, lighting designers would like to have the ability to bothvisualize and photometrically analyze these electric lighting anddaylighting system in a virtual environment as a design tool. Ideally,it should be possible (especially for theatrical and entertainmentlighting systems) to interactively control the lighting channel settingswhile interactively viewing the three-dimensional environment.

One prior art approach has been to specify a fixed virtual cameraposition and view direction, then generate a multiplicity of staticimages wherein for each image, exactly one of the lighting channels isset at full intensity while the other channels are disabled. Theseimages are then blended in real time during display to simulate theeffect of dimming the different lighting channels. This approach ismostly limited, however, to the representation of simple architecturaldimming system, as each lighting channel requires a full-resolutionstill image. A theatrical lighting system with hundreds to thousands oflighting channels would overwhelm the graphics capabilities of mostdesktop computers, not to mention the more limited capabilities oftablet computers and smartphones.

Another prior art approach relies on the real-time display capabilitiesof graphics hardware subsystems such as OpenGL from Khronos Group andDirect3D from Microsoft. Commercial software products such as, forexample, WYSIWIG from Cast Software, enable lighting designers tovisualize theatrical lighting systems with hundreds to thousands oflighting channels controlled by physical lighting consoles that areconnected to a desktop computer. The disadvantage of this approach isthat, while the 3D virtual environment can be interactively viewed, onlythe direct illumination from the luminaires can be modeled; indirectillumination due to interreflections from surfaces (which can oftenconstitute over half of the illumination in an enclosed environment) isnecessarily ignored. Further, the luminous intensity distributions ofthe luminaires—a critical aspect of professional lighting design—canonly be approximated using the real-time capabilities of the graphichardware subsystems. The resultant renderings are schematic rather thanphysically accurate.

There is therefore as need for a system and method of real-time dynamiclighting simulation, wherein a photometrically accurate representationof an environment can be displayed and analyzed in real time, while aplurality of hundreds to thousands of lighting channels may have theirintensity settings being continually varied and the user mayinteractively view the three-dimensional environment, all without theneed for ongoing global illumination calculations. The system and methodshould further be suitable for implementation on a variety of platforms,ranging from desktop computers to smartphones with bandwidth-limitedcommunications.

SUMMARY

A method performed by a lighting modelling system for displayingreal-time dynamic lighting simulations, the method comprising: encodinga texture map as a multiplicity of canonical radiosity solutions, eachrepresenting a lighting channel; storing the solutions in the texturememory of a graphics processing unit; generate a multiplicity oflighting channel intensity settings; store the lighting channelintensity settings; access the canonical radiosity solutions on aper-vertex basis with a vertex shader; multiply the associated vertexchannel colors by the lighting channel intensity settings; and sum theresultant colors to generate a vertex color for display.

The disclosed and/or claimed subject matter is not limited by thissummary as additional aspects are presented by the following writtendescription and associated drawings.

BRIEF DESCRIPTION OF DRAWINGS

FIG. 1A shows a block diagram of the predictive daylight harvestingsystem for buildings.

FIG. 1B shows a block diagram of the predictive daylight harvestingsystem for greenhouses.

FIG. 2 shows an example of a finite element representation of a combinedinterior and exterior environment.

FIG. 3 shows the transfer of luminous flux from exterior environmentelements through a window or opening to interior environment elements.

FIG. 4 shows a flowchart for operation of the predictive daylightharvesting system.

FIG. 5A shows a flowchart for reading the input data for buildings.

FIG. 5B shows a flowchart for reading the input data for greenhouses.

FIG. 6 shows a flowchart for calculating the optimal luminous andthermal environment.

FIG. 7A shows a flowchart for controlling the electric lighting,automated fenestration devices, and HVAC devices for buildings.

FIG. 7B shows a flowchart for controlling the electric lighting,automatic fenestration devices, HVAC devices, and supplemental carbondioxide distribution devices for greenhouses.

FIG. 8 shows a schematic representation of a predictive daylightingharvesting system.

FIG. 9 shows a block diagram of a predictive daylight harvesting system.

FIG. 10 is a flow diagram depicting an example method.

FIG. 11 is a schematic diagram depicting an example computing system.

FIG. 12 shows the projection of two virtual images with randomlyselected pixels onto a single depth buffer.

FIG. 13 shows the merger of two virtual images into a single virtualimage.

FIG. 14 shows the projection of a single virtual image through twopatches onto interior elements.

FIG. 15 is a flowchart depicting a method of transferring elementidentifiers through a transition surface.

FIG. 16 is a flowchart depicting another method of transferring elementidentifiers through a transition surface.

FIG. 17 is a flowchart depicting a method of geometric simplificationusing the marching cubes method.

FIG. 18 is a flowchart depicting a method of geometric simplificationusing Morton codes and adaptive cuts of k-d trees.

FIG. 19 illustrates the geometry of a glossy reflection.

FIG. 20 plots approximate and exact Blinn-Phong normalization termsversus the Blinn-Phong exponent.

FIG. 21 plots the approximate Blinn-Phong normalization term error.

FIG. 22 illustrates the geometry of translucent surface transmittance.

FIG. 23 shows example analytic bidirectional transmittance distributionfunctions.

FIG. 24 show the transfer of flux Φ from source patch S to receiverpatch R.

FIG. 25 show the process of form factor determination using Nusselt'sanalogy.

FIG. 26 shows the differential form factor dF between infinitesimalpatches dS and dR.

FIG. 27 projection of a planar surface onto the faces of a hemicube.

FIG. 28 shows the hemicube face cell geometry.

FIG. 29 shows the cubic tetrahedron geometry.

FIG. 30 shows the cubic tetrahedron face cell geometry.

FIG. 31 shows an example of daylight transmittance through a translucentwindow with a Blinn-Phong exponent of m=2.0.

FIG. 32 shows an example of daylight transmittance through a translucentwindow with a Blinn-Phong exponent of m=10.0.

FIG. 33 shows an example of daylight transmittance through a translucentwindow with a Blinn-Phong exponent of m=50.0.

FIG. 34 is a flowchart depicting the steps in diffusion properties of atransition surface material.

FIG. 35 shows the mean and first three principal component spectralreflectance distributions of 3,534 natural and synthetic materials.

FIG. 36 shows the approximation of four randomly-selected measuredspectral reflectance distributions by the weighted sums of six principalcomponents

FIG. 37 shows the approximation of four randomly-selected measuredspectral reflectance distributions by the weighted sums of threeprincipal components.

FIG. 38 shows the determination of a virtual photometer canonicalspectral irradiance distribution as the summation of direct spectralirradiance and indirect spectral irradiance.

FIG. 39 shows the multiplication of a light source complex spectralpower distribution by a virtual photometer canonical spectral irradiancedistribution to obtain a complex spectral irradiance distribution.

FIG. 40 is a flowchart depicting the determination of complex spectralirradiance distribution.

FIG. 41 shows the n canonical radiosity solutions of a virtualenvironment.

FIG. 42 shows a flowchart for the data storage of a multiplicity ofcanonical radiosity solutions in a two-dimensional texture map.

FIG. 43 shows a system for displaying virtual real-time dynamic lightingsystems using an external lighting console as an input source.

FIG. 44 shows a flowchart of a vertex shader program for displaying aninteractive dynamic lighting simulation in real time.

FIG. 45 shows a flowchart of a fragment shader program for displaying aninteractive dynamic lighting simulation in real time.

DETAILED DESCRIPTION

Electric lighting accounts for approximately 40 percent of all energyconsumed in modern buildings. Incorporating available daylight canreduce these annual energy costs by 40 to 50 percent using “daylightharvesting” techniques. The basic principle of daylight harvesting is tomonitor the amount of daylight entering an interior space and dim theelectric lighting as required to maintain a comfortable luminousenvironment for the occupants. Where required, motorized blinds andelectrochromic windows may also be employed to limit the amount ofdaylight entering the occupied spaces. Further energy savings can berealized through the use of occupancy sensors and personal lightingcontrols that operate in concert with the daylight harvesting system andare therefore considered integral thereto.

As a non-limiting example, the present invention provides a predictivedaylight harvesting system for buildings which can be implemented as amethod comprising the steps of: a) inputting data values regarding aplurality of variable building design parameters; b) calculating theeffects on a building's environmental characteristics based on the datavalues regarding a plurality of building design parameters; c) changingat least one of the data values regarding variable building designparameters; d) recalculating the effects on a building's environmentalcharacteristics based on the data values regarding a plurality ofbuilding design parameters. Buildings in the context of this patent areinclusive of any structure wherein the ambient environment may beaugmented through controlled means, including but not limited to light,heat, and humidity. As a result, buildings include, but are not limitedto, residential premises, commercial premises (including offices andwarehouses), aviculture facilities, agriculture facilities and animalhusbandry facilities. In addition to buildings in general, specificreference will be made to greenhouses and aquaculture facilities withoutlimiting the application of the predictive daylight harvesting systemdescribed herein.

Daylight harvesting is important wherever the amount of availabledaylight must be controlled or supplemented with electric lighting.Examples include hospitals and long-term care facilities where dailyexposure to natural lighting is beneficial, fish farms and aquaponicfacilities where overhead or underwater lighting is used to entrain thecircadian rhythms of fish, and free-range aviculture where birds requirea minimum amount of natural and electric lighting for optimum health andwell-being.

The parameters may be used to control any building environmental systemincluding lighting, heating, humidity or all of them together. Thisenables the selection of building design parameter settings to maximizeenergy savings while maintaining selected minimal occupant, includingplant or animal, requirements, for example, for heating, lighting,and/or humidity. Daylight harvesting is well suited to lighting controland to a lesser extent heating control, but the system is not limited tothat prime example of a use for the system.

Solar insolation is not necessarily a parameter to be used in thesystem, but it is a prime parameter to be considered in building design.In practice it would be preferred to include the steps of: a) measuringactual solar insolation and fine-tuning selected building designparameter settings to maximize energy savings while maintaining selectedminimal occupant, including plant or animal, requirements for heatingand lighting; b) analyzing a combination of historical weather data, insitu daylight measurements over time, and current weather predictions,and determining an optimal strategy for predictive daylight harvestingthat maximizes energy savings while maintaining selected minimaloccupant, including plant or animal, requirements for heating andlighting; c) analyzing occupant behavior, based on input from occupancyevent sensors and personal lighting control actions, or plant or animalgrowth and health based on input from plant sensors and/or manualgreenhouse operator input, and determining an optimal strategy fordaylight harvesting that maximizes energy savings while maintainingselected minimal occupant, including plant or animal, requirements forheating and lighting, based on predicted occupant behavior, or plant oranimal growth and health; and d) interacting with a building's HVACcontrol system and implementing an optimal strategy for maximizingenergy savings while maintaining selected minimal occupant requirements,including plant or animal, for heating, ventilation, air-conditioning,and lighting.

The variable building design parameters would include one or moredaylight photosensor locations, luminaire locations and illuminationlevels, occupancy sensor locations, temperature sensor locations andtemperatures, humidity sensor locations and humidity levels, glazingtransmittance characteristics, and thermal emissivity and thermal massof objects and surfaces in a building's interior environment for thepurpose of determining radiative heat transfer within the buildingenvironment due to solar insolation.

The variable greenhouse design parameters would include one or moredaylight photosensor locations, luminaire locations and illuminationlevels, temperature sensor locations and temperatures, humidity sensorlocations and humidity levels, soil moisture levels, carbon dioxide andoxygen gas concentrations, plant phytometrics, wind speed and direction,glazing transmittance characteristics, and thermal emissivity andthermal mass of objects and surfaces in a building's interiorenvironment for the purpose of determining radiative heat transferwithin the building environment due to solar insolation.

As another non-limiting example, in a basic implementation, thepredictive daylight harvesting method would have the following steps: a)receiving input data from at least one ambient condition sensor and atleast one information feed about anticipated solar conditions; b)calculating a luminous environment for a building based on the inputdata; c) generating messages based on output data about the calculatedluminous environment; and d) transmitting the messages via aninterconnect system to a building environmental control subsystem, inorder to maximize energy savings while maintaining selected minimaloccupant, including plant or animal, requirements for a building'senvironmental characteristics.

As yet another non-limiting example, the predictive daylight harvestingsystem can be implemented as an apparatus having a controller that: a)reads input data from a variety of sensors and information feeds, thesensors and feeds to include at least a plurality of sensors andinformation feeds from among the class of sensors and information feedsthat includes daylight photosensors, temperature sensors, occupancysensors, humidity sensors, soil moisture sensors, gas concentrationsensors, anemometers, timers, personal lighting controls, utility smartpower meters, weather report feeds, HVAC and energy storage controllers;b) calculates the effects of variable building design parameters onbuilding environment characteristics, such as light, heat, and humiditybalance and on energy management; and c) outputs building designparameter setting command signals, in order to maximize energy savingswhile maintaining selected minimal occupant, including plant or animal,requirements for the building environment characteristics respectively.The controller would read input data from a variety of sensors andinformation feeds, including but not limited to daylight photosensors,temperature sensors, occupancy sensors, humidity sensors, soil moisturesensors, gas concentration sensors, anemometers, timers, personallighting controls, utility power meters, weather report feeds, and otherenergy management systems, including HVAC and energy storagecontrollers. The controller would receive and process information aboutlight fixtures and light sources (luminaires) located in a building'sinterior environment, including photometric and electrical properties ofthe luminaires.

The controller may also read input data from and output command signalsand data to: a) variable power generation systems, including windturbines, solar panels, tidal power generators, run-of-river powergenerators, geothermal power generators, electrical batteries,compressed air and molten salt storage systems, biomass powergenerators, and other off-grid power sources; b) local heat sources,including geothermal and biomass fermenters; c) water recycling andreclamation systems, including rainfall collectors and cisterns,snowmelt, and atmospheric water generators; and d) air pollution controlsystems, including air filters and scrubbers.

The building control system would be enhanced by having the controllermaintain and access virtual representations of a building's exterior andinterior environments, including the geometry and material properties ofobjects that may significantly influence the distribution of daylightand artificial (for example, electrically-generated) luminous fluxwithin the environments, such as luminaires located in the interiorenvironment, including their photometric and electrical properties,daylight photosensors located in the interior and optionally exteriorenvironments, and occupancy sensors located in the interior environment.These virtual representations of building exterior and interiorenvironments would be accessed by the controller to perform calculationson the effects of solar insolation on building heat balance and energymanagement. Specifically, virtual representations of thermal emissivityand heat capacity (thermal mass) of objects and surfaces in a building'sinterior environment, would accessed by the controller for the purposeof determining radiative and conductive heat transfer within theenvironment due to solar insolation. Additionally, virtualrepresentations of occupants and their behaviors, including where theoccupants are likely to be located within a building's interiorenvironments at any given time and date, and the occupants' personallighting preferences, would be accessed by the controller for thepurpose of calculating optimal output settings for building designparameter setting command signals in order to maximize energy savingswhile maintaining selected minimal occupant requirements for heating andlighting.

Similarly, a greenhouse control system would be enhanced by having thecontroller maintain and access virtual representations of a greenhouse'sexterior and interior environments, including the geometry and materialproperties of objects that may significantly influence the distributionof daylight and artificial (for example, electrically-generated)luminous flux within the environments, such as luminaires located in theinterior environment, including their photometric and electricalproperties, daylight photosensors, and humidity sensors located in theinterior and optionally exterior environments, soil moisture and gasconcentration sensors located in the interior environment, and windspeed and direction sensors located in the exterior environment. Thesevirtual representations of building exterior and interior environmentswould be accessed by the controller to perform calculations on theeffects of solar insolation on greenhouse heat balance and energymanagement. Specifically, virtual representations of thermal emissivityand heat capacity (thermal mass) of objects and surfaces in agreenhouse's interior environment, would accessed by the controller forthe purpose of determining radiative and conductive heat transfer withinthe environment due to solar insolation. Additionally, virtualrepresentations of plant growth and health, including the plants'photosynthetic and photomorphogenetic processes, temperature, humidityand other environmental requirements, would be accessed by thecontroller for the purpose of calculating optimal output settings forgreenhouse design parameter setting command signals in order to maximizeenergy savings while maintaining selected minimal plant requirements forheating and lighting.

Optimally, the building control system would include a fuzzy-logicinference engine that learns weather patterns and occupant usagepatterns and preferences, and the controller would maintain amathematical model of sky conditions, historical weather data, and adatabase of weather patterns and occupant usage patterns andpreferences, continually reads data from external input andcommunication devices, calculates the optimal settings for luminairesand fenestration devices, and controls luminaires and fenestrationdevices to achieve maximal energy savings while providing an interiorluminous and thermal environment that is consistent with predefinedgoals and occupant preferences.

Similarly, an optimal greenhouse control system would include afuzzy-logic inference engine that learns weather patterns and optimalplant and animal growth and health parameters, and the controller wouldmaintain a mathematical model of sky conditions, historical weatherdata, and a database of weather patterns and plant short-term andlong-term environmental requirements, continually reads data fromexternal input and communication devices, calculates the optimalsettings for luminaires and fenestration devices, and controlsluminaires, fenestration devices, humidity control devices, andsupplemental carbon dioxide distribution devices to achieve maximalenergy savings while providing an interior luminous and thermalenvironment that is consistent with predefined goals and plantrequirements.

In one elementary form, the predictive daylight harvesting system wouldalso comprise: a) at least one controller that reads input data from avariety of sensors and information feeds, and that includes anartificial intelligence engine; b) at least one ambient condition sensorand at least one information feed; and c) an interconnect systemoperatively coupling the controller to the sensor and the informationfeed, configured to provide output data suitable for dimming orswitching luminaires and operating automated fenestration and otherenvironmental devices.

The controller may further maintain communication with other buildingautomation subsystems, including but not limited to HVAC and energystorage systems. It may also maintain communication with externalsystems such as electrical power utilities and smart power grids.

In a preferred mode of operation, the controller would continually readdata from the external input and communication devices, calculate theoptimal settings for the luminaires, fenestration, and otherenvironmental control devices, and control those devices to achievemaximal annual energy savings while providing an interior luminous andthermal environment that is consistent with predefined goals andoccupant preferences or plant or animal requirements. The “what-if”scenarios capability of the invention deriving from its simulating avirtual building interior environment on a regular basis (for example,hourly) enable a physical daylight harvesting controller system to bedesigned (for example, including an optimal layout of daylightphotosensors) and programmed accordingly. The controller can thenfurther access the virtual representation during operation to refine itsbehavior in response to the building performance by means of “what-if”simulations.

The present invention is herein described more fully with reference tothe accompanying drawings, in which embodiments of the invention areshown. This invention may, however, be embodied in many different forms,and should not be construed as limited to the embodiments describedherein. Rather, these embodiments are provided so that this disclosurewill be thorough and complete, and will fully convey the scope of theinvention to those skilled in the art. As used in the presentdisclosure, the term “diffuse horizontal irradiance” refers to that partof solar radiation which reaches the Earth as a result of beingscattered by the air molecules, aerosol particles, cloud particles, orother airborne particles, and which is incident upon an unobstructedhorizontal surface, the term “direct normal irradiance” refers to solarradiation received directly from the Sun on a plane perpendicular to thesolar rays at the Earth's surface, the term “solar insolation” refers tothe solar radiation on a given surface, and the terms “illuminance,”“irradiance,” “luminous exitance,” “luminous flux,” “luminousintensity,” “luminance,” “spectral radiant flux,” “spectral radiantexitance,” and “spectral radiance” are commonly known to those skilledin the art.

As used in the present disclosure, the term “photosynthetically activeradiation” (abbreviated “PAR”) refers to the spectral range ofelectromagnetic radiation (either solar or electrically-generated) from400 nm to 700 nm that photosynthetic organisms are able to use in theprocess of photosynthesis. It is commonly expressed in units ofmicromoles per second (μmol/sec).

As used in the present disclosure, the term “photosynthetic photon fluxdensity” (abbreviated “PPFD”) refers to the density of photosyntheticphotons incident upon a physical or imaginary surface. It is commonlymeasured in units of micromoles per square meter per second(μmol/m2-sec) with a spectrally-weighted “quantum” photosensor.

As used in the present disclosure, the term “ultraviolet radiation”refers to the spectral range of electromagnetic radiation (either solaror electrically-generated) from 280 nm to 400 nm. It may be furtherreferred to as “ultraviolet A” (abbreviated “UV-A”) with spectral rangeof 280 nm to 315 nm, and “ultraviolet B” (abbreviated “UV-B”) with aspectral range of 315 nm to 400 nm. It may be expressed in units ofmicromoles per second (μmol/sec) or watts.

As used in the present disclosure, the term “visible light” refers tothe spectral range of electromagnetic radiation (either solar orelectrically-generated) from 400 nm to 700 nm.

As used in the present disclosure, the term “infrared radiation” refersto the spectral range of electromagnetic radiation (either solar orelectrically-generated) from 700 nm to 850 nm for horticulturalpurposes, and from 700 nm to 2800 nm for solar insolation purposes. Itmay be expressed in units of micromoles per second (μmol/sec) or watts.

As used in the present disclosure, the term “environment” may refer to avirtual environment comprised of a finite element model or similarcomputer representation with virtual sensors and control devices, or aphysical environment with physical sensors and control devices.

As used in the present disclosure, the term “luminous environment”refers to a physical or virtual environment wherein visible light,ultraviolet radiation, and/or infrared radiation, is distributed acrossmaterial surfaces.

As used in the present disclosure, the term “solar heat gain” refers tothe increase in temperature in an object that results from solarradiation.

Electric lighting accounts for approximately 40 percent of all energyconsumed in modern buildings. It has long been recognized thatincorporating available daylight can reduce these annual energy costs by40 to 50 percent using “daylight harvesting” techniques.

Daylight harvesting is also an essential feature of agriculturalgreenhouses, whose primary function is to provide a controlledenvironment for growing vegetables or flowers. Optimizing daylight usageminimizes supplemental lighting and heating costs while maintainingoptimal growing conditions.

The basic principle of daylight harvesting is to monitor the amount ofdaylight entering an interior space and dim the electric lighting asrequired to maintain a comfortable luminous environment for theoccupants. Where required, motorized blinds and electrochromic windowsmay also be employed to limit the amount of daylight entering theoccupied spaces.

The same principle applies to greenhouses, where it is necessary toprovide supplemental electric lighting and heating or ventilation asrequired to maintain optimal growing conditions. Where required,motorized shading devices and electrochromic windows may also beemployed to limit the amount of daylight entering the greenhouse spaces.

Further energy savings can be realized through the use of occupancysensors and personal lighting controls that operate in concert with thedaylight harvesting system and are therefore considered integralthereto. At the same time, the building occupants' comfort andproductivity must be taken into consideration by for example limitingvisual glare due to direct sunlight.

Yet further energy savings can be realized through the use oftransparent or semi-transparent solar panels, particularly those whichabsorb and convert to electricity ultraviolet and/or infrared radiationwhile remaining substantially transparent to visible light. Theperformance of such panels when used for building glazing therebybecomes integral to the performance of the daylight harvesting system.

Daylight harvesting is also related to solar insolation management inthat the infrared solar irradiation entering interior spaces is absorbedby the floors, walls and furniture as thermal energy. This influencesthe building heat and humidity balance that must be considered in thedesign and operation of heating, ventilation and air conditioning (HVAC)systems. For example, the energy savings provided by electric lightdimming may need to be balanced against the increased costs of HVACsystem operation needed to maintain building occupant's comfort andproductivity or plant growth and health. The use of temperature sensorsand humidity sensors are therefore also considered integral to adaylight harvesting system.

For greenhouses, the amount of daylight entering the greenhouse space isdirectly related to plant photosynthesis and water evaporation. Furtherenergy savings and optimization of plant growth and health can thereforebe realized through the use of temperature sensors, humidity sensors,soil moisture sensors, carbon dioxide and oxygen concentration sensors,and sensors for directly monitoring plant growth and health (forexample, U.S. Pat. No. 7,660,698 and U.S. Pat. No. 8,836,504). Suchsensors operate in concert with the daylight harvesting system and aretherefore considered integral thereto.

Related to solar insolation management is the design and operation ofsolar energy storage systems for high performance “green” buildings,wherein thermal energy is accumulated by means of solar collectors andstored in insulated tanks or geothermal systems for later use in spaceheating. Similarly, electrical energy may be accumulated by means ofsolar photovoltaic panels or wind power plants and stored in batteriesfor later use in operating electric lighting. The decision of whether tostore the thermal and/or electrical energy or use it immediately isdependent in part on daylight harvesting through the operation ofmotorized blinds and other active shading devices that affect solarinsolation of interior surfaces.

The operation of such systems may therefore benefit from weather dataand other information that is shared by the daylight harvesting system.More generally, the design and operation of a daylight harvesting systemis advantageously considered a component of an overall buildingautomation system that is responsible for energy usage and conservation.

Annual lighting energy savings of 40 to 50 percent are possible withwell-designed electric lighting systems incorporating daylightharvesting, even as standalone controllers that function independentlyof the HVAC systems. However, it has also been shown that roughly halfof all installed daylight harvesting systems do not work as designed.These systems do not provide significant energy savings, and in manycases have been physically disabled by the building owners due tounsatisfactory performance.

The underlying problem is that the performance of these systems isdetermined by many parameters, including lamp dimming (continuous versusbi-level switched), photosensor placement and orientation, photosensorspatial and spectral sensitivity, luminaire zones, timer schedules,occupancy sensors, interior finishes, window transmittances, exteriorand interior light shelves, motorized blinds and other shading devices,and personal lighting controls. In addition, the presence of exterioroccluding objects such as other buildings, large trees and surroundinggeography need to be considered, as do the hour-by-hour weatherconditions for the building site.

A similar underlying problem exists with greenhouses in that theperformance of greenhouse control systems is determined by manyparameters, including lamp switching, photo sensor placement andorientation, photosensor spatial and spectral sensitivity, luminairezones, timer schedules, humidity sensors, soil moisture sensors, gassensors, and plant monitoring sensors, interior finishes, exterior andinterior light shelves and similar reflectors, motorized shadingdevices, heating and ventilation devices, humidity control devices,supplemental carbon dioxide distribution devices, and greenhouseoperator controls. In addition, the presence of exterior occludingobjects such as other buildings, large trees and surrounding geographyneed to be considered, as do the hour-by-hour weather conditions for thegreenhouse building site. Wind direction and speed is also a concern,due to the typically high rate of conductive heat loss through thegreenhouse glazing system.

Given all this, there are at present no suitable design tools for thearchitect or lighting designer to simulate the performance of suchcomplex systems. In particular, the current state-of-the-art in lightingdesign and analysis software requires ten of minutes to hours ofcomputer time to simulate a single office space for a given time, date,and geographic location. Architects and lighting designers have nochoice but to calculate the lighting conditions at noon on the winterand summer solstices and spring and fall equinoxes, then attempt toestimate the performance of their daylight harvesting designs for theentire year based solely on this limited and sparse information.

Existing daylight harvesting systems are therefore designed according tobasic rules of thumb, such as “locate the photosensor away from directsunlight” and “dim the two rows of luminaires closest to thesouth-facing windows.” There is no means or opportunity to quantify howthe system will perform when installed.

There are software tools for building HVAC systems design that takesolar insolation into account, but they do so in the absence of designinformation concerning daylight harvesting and electric lighting. Themechanical engineer typically has to assume so many watts per unit areaof electric lighting and design to worst-case conditions, with noconsideration for energy savings due to daylight harvesting.

Once the daylight harvesting system has been installed, it must becommissioned. This typically involves a technician visiting the buildingsite once on a clear day and once at night to calibrate the photosensorresponses. This is mostly to ensure that the system is working, as thereis no means or opportunity to adjust the feedback loop parametersbetween the photosensors and the dimmable luminaires for optimumperformance. Indeed, many existing daylight harvesting systems operatein open loop mode without any feedback, mostly for the sake of designsimplicity.

There is therefore a need for a method and apparatus to accomplish fiverelated tasks: First, there is a need for a method whereby an architector lighting designer can fully and interactively simulate theperformance of a daylight harvesting design. That is, the method shouldaccount for multiple system parameters that may influence the systemperformance, including the effects of solar insolation on building heatbalance and energy management. The method should simulate the overallsystem such that an hour-by-hour simulation for an entire year can begenerated in at most a few minutes. Given this interactive approach, anarchitect or lighting designer can experiment with many differentdesigns or variations on a design in order to determine which designmaximizes annual building energy savings while respecting buildingoccupants' comfort or plant or animal growth and health requirements.

Second, there is a need for an apparatus that uses the aforesaid methodto physically implement the simulated system, and which has the abilityto self-tune its initial parameters (as determined by the simulation) inorder to maximize annual energy savings.

Third, there is a need for an apparatus that can analyze a combinationof historical weather data, in situ daylight measurements over time,current weather predictions, and other information using the aforesaidmethod to autonomously determine an optimal strategy for predictivedaylight harvesting that maximizes annual energy savings.

Fourth, there is a need for an apparatus that can analyze occupantbehavior, including occupancy sensor events and personal lightingcontrols actions such as luminaire dimming and switching and setting ofpersonal time schedules, in buildings to autonomously determine anoptimal strategy for daylight harvesting that maximizes annual energysavings based on predicted occupant behavior.

Similarly, there is a need for an apparatus that can analyze plantgrowth and health as determined by automatic plant growth sensors ormanual greenhouse operator input, and autonomously determine an optimalstrategy for daylight harvesting that maximizes annual energy savingsbased on predicted plant growth and health.

Fifth, there is a need for a daylight harvesting apparatus that cancoordinate its operation with HVAC systems, solar energy storagesystems, and other high performance “green” building managementtechnologies that are influenced by solar insolation for the purposes ofmaximizing annual energy savings and building occupants' comfort andproductivity or plant growth and health.

FIG. 1A shows an apparatus for enabling a predictive daylight harvestingsystem for a building. As shown, it is logically, but not necessarilyphysically, comprised of three components: 1) inputs 10 to 40; 2)controller 45, user interface 50 and communications port 55; and 3)outputs 60 and 65.

FIG. 1B shows an apparatus for enabling a predictive daylight harvestingsystem for a greenhouse. As shown, it is logically, but not necessarilyphysically, comprised of three components: 1) inputs 100 to 140; 2)controller 145, user interface 150 and communications port 155; and 3)outputs 160 to 175.

Inputs

Referring to FIG. 1A, the inputs include daylight photosensors 10,occupancy or vacancy sensors 15, timers 20, personal lighting controls25, utility power meters 30, weather report feeds 35, and optionaltemperature sensors 40.

Referring to FIG. 1B, the inputs include daylight photosensors 100,temperature sensors 105, timers 110, humidity sensors 115, and optionalwind speed and direction sensors 120, soil moisture sensors 125, carbondioxide and/or oxygen sensors 130, utility smart power meters 135, andweather report feeds 140.

For the purposes of this application, a photosensor 10 or 100 includesany electronic device that detects the presence of visible light,infrared radiation (IR), and/or ultraviolet (UV) radiation, includingbut not limited to photoconductive cells and phototransistors. Thephotosensor may be spectrally weighted to measure units of illuminance,irradiance, PPFD, ultraviolet power flux density, infrared power fluxdensity, or spectral radiant flux.

For the purposes of this application, a gas concentration sensor 130includes any electronic device that measures the partial pressure of agas, including but not limited to, carbon dioxide and oxygen.

Spectral radiant flux may be measured with a spectroradiometer. Suchinstruments are useful for greenhouse application in that whileknowledge of photosynthetically active radiation (PAR) and PPFD isuseful for predicting plant growth and health on a per-species basis, itonly measures visible light integrated across the spectral range of 400nm to 700 nm. It is known however that ultraviolet radiation may inducechanges in leaf and plant morphology through photomorphogenesis, thatinfrared radiation influences seed germination, flower induction, andleaf expansion, and that the presence of green light is important tomany plants. The response of different plant species to spectral powerdistributions becomes important as multi-color LED luminaires replacehigh-pressure sodium and other lamp types with fixed spectral powerdistributions.

Similarly, for the purposes of this application, an occupancy or vacancysensor 15 includes any electronic device that is capable of detectingmovement of people, including but not limited to pyroelectric infraredsensors, ultrasonic sensors, video cameras, RFID tags, and securityaccess cards.

Commercial daylight photosensors 10 or 100 may feature user-selectablesensitivity ranges that may be specified in foot-candles or lux(illuminance), watts/m2 (irradiance), or μmol/m2-sec (PPFD).Additionally, the photosensor may include a field of view which may bespecified in degrees. For the purposes of the present invention thespatial distribution of the photosensor sensitivity within its field ofview, the spectral responsivity of the photosensor in the visible,infrared and ultraviolet portion of the electromagnetic spectrum, andthe transfer function of the electrical output from the photosensor aregermane. These device characteristics enable the present invention tomore accurately simulate their performance, though the exclusion of oneor more of these parameters will not prevent the system fromfunctioning.

Daylight photosensors may be positioned within an interior environment(such as the ceiling of an open office) such that they measure theaverage luminance of objects (such as floors and desktops) within theirfield of view. Positioning of the photosensors may be to locate thephotosensor away from direct sunlight. However, so-called “dual-loop”photosensors as disclosed in U.S. Pat. No. 7,781,713 and U.S. Pat. No.7,683,301 may be positioned in for example skywells and skylights tomeasure both the average luminance of objects below and the directsunlight and diffuse daylight incident upon the photosensors from above.

An alternative is to locate a daylight photosensor in each luminaire. Inthis example, the luminaire may be dimmed according to how much ambientlight is detected by its photosensor.

Commercial occupancy and vacancy sensors 15 may employ pyroelectric,ultrasonic, or optical detection technologies, or a combination thereof.It can be difficult to accurately characterize the performance of thesedevices in enclosed spaces, as they may be influenced by the thermalemissivity of building materials for reflected far-infrared radiationand the acoustic reflectance of building materials for reflectedultrasonic radiation. In the present invention, occupancy sensors workbest according to line-of-sight operation within their environments andin accordance with their sensitivity in terms of detection distance.Other motion detection techniques may also be employed.

Timers 20 or 110 may be implemented as external devices that areelectrically or wirelessly connected to the building controller 45 orgreenhouse controller 145, or they may be implemented within thehardware and software of said controller and accessed through thebuilding controller's user interface 50 or greenhouse controller's userinterface 150.

Personal lighting controls 25 may be implemented as for example handheldinfrared or wireless remote controls or software programs executed on adesktop computer, a laptop computer, a smartphone, a personal digitalassistant, or other computing device with a user interface that can beelectrically or wirelessly connected to the building controller 45,including an Internet connection, on a permanent or temporary basis. Thecontrols optionally enable the occupant or user to specify for examplepreferred illuminance levels in a specific area of the interiorenvironment (such as, for example, an open office cubicle or privateoffice), to control selected motorized blinds or electrochromic windows,to influence (such as, for example, by voting) the illuminance levels ina shared or common area of the interior environment, to specify minimumand maximum preferred illuminance levels, to specify the time rate ofillumination level increase and decrease (“ramp rate” and “fade rate”),to specify time delays for occupancy sensor responses, and to specifyindividual time schedules.

Utility smart power meters 30 or 135 can provide real-time informationon power consumption by buildings, in addition to information onvariable power consumption rates that may change depending on theutility system load and policies. Such meters may be electrically orwirelessly connected to the building controller 45 or greenhousecontroller 145, including an Internet connection.

Real-time weather report feeds 35 or 140 are widely available includingon the Internet. These feeds can be connected to the building controller45 or greenhouse controller 145 via a suitable electrical or wirelessconnection.

Geographic Information Systems (GIS) data available through varioussources can additionally provide relevant environmental data, againusing a feed connected to the build controller 45 or greenhousecontroller 145 via a suitable electrical or wireless connection.

Temperature sensors 45 or 105 may be employed to measure room or spacetemperatures if such information is not available from an external HVACcontroller or building automation system (not shown) in communicationwith building controller 45 or greenhouse controller 145.

Controller

In an embodiment, the building controller 45 or greenhouse controller145 is a standalone hardware device, advantageously manufactured to thestandards of standalone industrial computers to ensure continuous andreliable operation. It can however be implemented as a module of alarger building automation system.

In another embodiment, the building controller 45 or greenhousecontroller 145 is comprised of a standalone hardware device, such as forexample a communications hub, that is in communication with an externalcomputing device, such as for example a centralized server, a networkedremote computer system, or a cloud-based software-as-a-service.

The building controller may further comprise a user interface 50 and oneor more communication ports 55 for communication with operators andexternal systems, including HVAC controllers, energy storage systemcontrollers, building automation systems, and geographically remotedevices and systems (not shown). Similarly, the greenhouse controller145 may further comprise a user interface 150 and one or morecommunication ports 155 for communication with operators and externalsystems, including HVAC controllers, energy storage system controllers,building automation systems, and geographically remote devices andsystems (not shown).

The operation of the building controller 45 or greenhouse controller 145is as disclosed below following a description of the outputs.

Outputs

Referring to FIG. 1A, building controller 45 provides electrical signalsto the dimmable or switchable luminaires 60 and optionally automatedfenestration devices 65, such as for example motorized window blinds andelectrochromic windows whose transmittance can be electricallycontrolled. Said electrical signals may include analog signals, such asfor example the industry-standard 0-10 volt DC or 4-20 milliamp signals,and digital signals using a variety of proprietary and industry-standardprotocols, such as DMX512 and DALI. The connections may be hard-wiredusing for example an RS-485 or derivative connection, an Ethernetconnection, or a wireless connection, such as for example Bluetooth,Zigbee, 6LoWPAN, or EnOcean.

Referring to FIG. 1B, greenhouse controller 145 provides electricalsignals to the dimmable or switchable luminaires 160, and optionallyautomated fenestration devices 165, such as for example motorized windowblinds and electrochromic windows whose transmittance can beelectrically controlled, in addition to optional humidity controldevices 170 and carbon dioxide distribution devices 175. Said electricalsignals may include analog signals, such as for example theindustry-standard 0-10 volt DC or 4-20 milliamp signals, and digitalsignals using a variety of proprietary and industry-standard protocols,such as DMX512 and DALI. The connections may be hard-wired using forexample an RS-485 or derivative connection, an Ethernet connection, or awireless connection, such as for example Bluetooth, Zigbee, 6LoWPAN, orEnOcean.

Controller Operation

In one embodiment, the building controller 45 or greenhouse controller145 maintains or accesses a three-dimensional finite elementrepresentation of the exterior and interior environments, such as isshown for example in FIG. 2. The representation is comprised of a set ofgeometric surface elements with assigned surface properties such asreflectance, transmittance, and color, and one or more electric lightsources (“luminaires”) with assigned photometric data, such as luminousintensity distribution. For the purposes of solar insolation analysis,thermal emissivity and heat capacity properties may also be included.

Radiative flux transfer or ray tracing techniques can be used to predictthe distribution of luminous flux or spectral radiant flux emitted bythe luminaires within the interior environment due to interreflectionsbetween surface elements, as disclosed in for example U.S. Pat. No.4,928,250. An example of a commercial lighting design and analysissoftware products that employs such techniques is AGi32 as manufacturedby Lighting Analysts Inc. (Littleton, Colo.). As disclosed by U.S. Pat.No. 4,928,250, said techniques apply to both visible and invisibleradiation, including the distribution of infrared radiation due to solarinsolation.

For daylighting analysis, the representation of the exterior environmentmay include other buildings and topographical features that may occludedirect sunlight and diffuse daylight from entering the windows andopenings of the interior environment. Similarly, said buildings andtopographical features may reflect direct sunlight and diffuse daylightinto the interior environment via the windows and openings.

In a rural setting, for example, the site topography may be extractedfrom a geographic information system and represented as a set of finiteelements. In an urban setting, the geometry and material properties ofnearby buildings may be extracted from a series of photographs, such asare available for example from Google Street View.

Exterior Environments

For exterior environments, the light sources are direct solar irradianceand diffuse irradiance from the sky dome. Given the direct normal anddiffuse horizontal irradiance measurements for a given time of day andJulian date, and the geographical position of the site in terms oflongitude and latitude, the Perez Sky or similar mathematical model maybe used to accurately estimate the distribution of sky luminance for anyaltitude and azimuth angles. The irradiance measurements may be made insitu or derived from historical weather data such as from the TypicalMeteorological Year (TMY) database for the nearest geographic location.(Infrared direct normal and diffuse horizontal irradiance measurementsare also often available.)

The Perez sky model does not include spectral power distributions fordaylight. However, a validated extension to the Perez sky model thatincludes the ability to model ultraviolet radiation can be used topredict daylight spectral power distributions.

The Perez sky model also does not include the effects of low-levelatmospheric phenomena, including naturally-occurring fog, dust, andman-made smog, on direct solar irradiance and diffuse irradiance fromthe sky dome and the spectral power distribution of daylight. However,additional mathematical models, historical weather data, and in situmeasurements can be used to account for these effects.

Typical Meteorological Year (TMY3) database records also provide hourlyreadings for cloud coverage, dry-bulb temperature, dew-pointtemperature, relative humidity, wind direction, and wind speed, andother information that can be incorporated into the simulation modelingfor greenhouse controllers.

Using radiative flux transfer techniques, the distribution of luminousor spectral radiant flux due to direct sunlight and diffuse daylightwithin the exterior environment can be predicted. In an embodiment, thesky dome is divided into a multiplicity of horizontal bands with equalaltitudinal increments. Each band is then divided into rectangular “skypatches” at regular azimuthal increments such that each patch hasroughly equal area. An example of such a division that yields 145patches is commonly referred to as the Tregenza sky dome division. Apreferred division however may be obtained by recursive subdivision ofan octahedron to yield 256 approximately equal-area patches.

In the same embodiment, the daily solar path across the sky domethroughout the year is similarly divided into a multiplicity of bandsthat parallel the solar path for selected dates, such as for exampleDecember 21, February 01/November 09, February 21/October 21, March09/October 5, March 21/September 21, April 6/September 7, April21/August 21, May 11/August 04, and June 21. The multiplicity of bandswill span 47 degrees, regardless of the site latitude. Each band is thendivided into rectangular “solar patches” at regular time intervals (forexample, one hour) such that each patch has roughly equal area. Asubdivision with for example nine bands and one-hour intervals will becomprised of 120 solar patches.

In one embodiment, spectral radiant flux is represented as three or moreseparate color bands that are identified for example as “red,” “green,”and “blue” in correspondence with their perceived colors, and iscommonly practiced in the field of computer graphics. In anotherembodiment, more color bands may be employed. For example, the spectralresponsivity of daylight photosensors may extend from ultraviolet toinfrared wavelengths. It may therefore be advantageous to representspectral radiant flux as for example 50 nanometer-wide color bands from300 nm to 1200 nm. In yet another embodiment, a single radiant fluxvalue may suffice, such as for example to represent infrared radiationfor solar insolation analysis.

For greater accuracy in representing the sky dome luminancedistribution, a sky dome division may alternately be chosen such thatthe differences between adjacent patch luminances are minimized. Forexample, the sun will traverse a 47-degree wide band of the sky domeover a period of one year. Smaller patches within this band may beemployed to more accurately represent the sky dome luminance for anygiven time and Julian date.

According to prior art as implemented for example by AGi32, each skydome patch represents an infinite-distance light source that illuminatesthe exterior environment with a parallel light beam whose altitudinaland azimuthal angles are determined by the patch center, and whoseluminous intensity is determined by the Perez sky model. (Other skymodels, such as the IES and CIE Standard General Sky, may also beemployed to predict the sky dome's spatial luminance distribution.)Similarly, each solar patch represents an infinite-distance light sourcethat illuminates the exterior environment with a parallel light beamwhose altitudinal and azimuthal angles are determined by the patchcenter, and whose luminous intensity is determined by the Perez skymodel. Once the luminous flux contributions from the direct sunlight anddiffuse daylight have been calculated and summed, radiative fluxtransfer techniques such as those disclosed in U.S. Pat. No. 4,928,250can be employed to calculate the distribution of luminous or spectralradiant flux due to interreflections between surface elements in theexterior environment.

For greenhouse applications, it is also possible to model crop canopiesfor use with radiative flux transfer techniques such as those disclosedin U.S. Pat. No. 4,928,250. Such models can be updated as the greenhousecrops grow and mature.

As will be known to those skilled in the arts of thermal engineering orcomputer graphics, each surface element in the finite elementrepresentation is assigned a parameter representing the luminous orspectral radiant flux that it has received but not yet reflected and/ortransmitted (the “unsent flux”), and another parameter representing itsluminous or spectral radiant exitance. (Infrared and ultraviolet radiantflux and exitance may also be considered without loss of generality.) Ateach iteration of the radiative flux transfer process, the unsent fluxfrom a selected element is transferred to all other elements visible tothat element. Depending on the reflectance and transmittance propertiesof each element, some of the flux it receives is reflected and/ortransmitted; the remainder is absorbed. The flux that is not absorbed isadded to both its unsent flux and luminous exitance parameters. Thetotal amount of unsent flux thus decreases at each iteration, and so theradiative flux transfer process converges to a “radiosity” solution,wherein the luminous exitance or spectral radiant exitance of everysurface element is known.

Canonical Radiosity

In a novel contribution of the present invention, this approach isextended by assigning a multiplicity of ‘n’ unsent flux parameters and‘n’ luminous or spectral radiant exitance parameters to each exteriorenvironment element, where ‘n’ is the total number of sky and solarpatches. Each sky patch and solar patch ‘i’ is assigned unit luminanceor spectral radiance, and its contribution to the luminous or spectralradiant flux of each exterior element is saved in its ‘i’th unsent fluxand luminous or spectral radiant exitance parameters.

Once the luminous or spectral radiant flux contributions from thediffuse daylight (but not direct sunlight) have been calculated,radiative flux transfer techniques can be employed to calculate thedistribution of luminous or spectral radiant flux due tointerreflections between surface elements in the exterior environmentfor each sky and solar patch. The result is the generation of ‘n’separate radiosity solutions for the exterior environment. Because thesky and solar patch luminances were not considered, these are referredto as ‘canonical’ radiosity solutions.

A particular advantage of this approach is that approximately 95 percentof the computational time needed to calculate the distribution ofluminous flux due to interreflections between surface elements isdevoted to calculating the “form factors” (i.e., geometric visibility)between elements, as disclosed in U.S. Pat. No. 4,928,250. Thus, if theprior art approach requires ‘x’ amount of time to calculate a singleradiosity solution, the present approach requires ‘x’+‘y’ time, where‘y’ is typically less than five percent. (The remainder of thecomputational time is consumed by other “housekeeping” activitiesrelated to manipulating the geometric and material data.)

Once the ‘n’ canonical radiosity solutions have been calculated, theradiosity solution for any sky luminance distribution for a specifiedtime and date can be calculated as a weighted sum of the canonicalradiosity solutions, where the weights are the luminances of theassociated sky and solar patches as predicted by the chosen sky model.Such solutions may be calculated in milliseconds, as opposed to minutesto hours for prior art approaches.

The present invention employs an approach wherein the solar altitudinaland azimuthal angles are used to determine the altitudinal and azimuthalangles of the four closest solar patch centers. The radiosity solutionfor the direct sunlight is then calculated as the weighted sum of thefour canonical radiosity solutions for these patches, with the weightsbeing the luminous intensities of the direct sunlight at the respectivealtitudinal and azimuthal angles.

Geometry Simplification

A disadvantage to assigning one or a multiplicity of ‘n’ unsent fluxparameters and ‘n’ luminous or spectral radiant exitance parameters toeach environment element, where ‘n’ is the total number of sky and solarpatches, is that complex environments with tens or hundreds of thousandsof polygonal elements may require excessive amount of memory. To addressthis issue, selected objects comprised of many elements (such as, for,example, furniture in an office) may be represented by geometricallysimplified objects with fewer elements for the purposes of canonicalradiosity calculations.

In a preferred embodiment, the elements of a selected object are firstinserted into an octree or similar spatial subdivision data structure.The intersection of the object surface elements with the octree nodeedges are then determined.

A simplified geometric representation consisting of possibly disjointtriangles is then constructed from the octree data structure using thewell-known “marching cubes” algorithm. There are many variations on thealgorithm that attempt to address deficiencies in the original algorithmdisclosed in U.S. Pat. No. 4,710,876, “System and Method for the Displayof Surface Structures Contained within the Interior Region of a SolidBody”.

A disadvantage of the marching cubes algorithm is that while it reducesthe geometric complexity of objects, it still tends to generate manysmall triangular elements on the scale of the dimensions of the octreenodes. It is therefore necessary to perform further geometricsimplification of the resultant object by merging nearly coplanartriangular elements into larger polygonal elements.

In one embodiment, after the canonical radiosity calculations have beenperformed using the simplified geometric objects, it is necessary todetermine the radiosity solution for any sky luminance distribution fora specified time and date that includes the original object geometry.This can be accomplished by using various well-known “final gather”techniques commonly used in global illumination methods.

Referring to FIG. 17, the method comprises: inserting the polygonalelements of a geometric object into an octree or similar spatialsubdivision data structure (Step 1710); determining the intersection ofthe polygonal elements with the octree node edges (Step 1720);constructing a simplified representation of the geometric object usingthe marching cubes algorithm (Step 1730); further simplifying therepresentation using the coplanar polygon merge algorithm (Step 1740);calculating a canonical radiosity solution using the geometricallysimplified objects (Step 1750); and determining a radiosity solution forany sky luminance distribution for a specified time and date thatincludes the original object geometry using a final gather technique(Step 1760).

In another embodiment, the vertices of the elements of a selected objectare sorted according to their Morton codes to define a set of leaf nodeswith error quadrics. (As will be known to those skilled in the art,Morton codes map three-dimensional vertex coordinates into aone-dimensional array whose ordering is equivalent to the depth-firsttraversal of an octree.) The cumulative sum of these Morton codes,termed a Morton integral, is then used to construct a k-d tree that isequivalent to a bounding volume hierarchy spatial subdivision datastructure. This tree is then traversed using the error quadrics for eachinternal node to obtain an adaptive cut that represents the vertices ofa simplified geometric object representative of the selected object,with the vertices remeshed to form the object surfaces.

Referring to FIG. 18, the method comprises: calculating the vertex errorquadrics for each object (Step 1810); sorting the vertices of eachelement according to their Morton codes (Step 1820); calculating thecumulative sum of the Morton codes for each object to generate a Mortonintegral (Step 1830); constructing a k-d tree of the Morton codes andtheir associated vertices (Step 1840); determining an adaptive cut ofthe k-d tree based on the vertex error quadrics as a representation ofthe vertices of the simplified geometric object (Step 1850); remeshingthe vertices to form object surfaces (Step 1860); further simplifyingthe representation using the coplanar polygon merge algorithm (Step1870); calculating a canonical radiosity solution using thegeometrically simplified objects (Step 1880); and determining aradiosity solution for any sky luminance distribution for a specifiedtime and date that includes the original object geometry using a finalgather technique (Step 1890).

Virtual Photosensors

A disadvantage of the radiative flux transfer techniques is that whilethe radiosity solution may provide the luminance and approximatespectral radiant exitance at any selected point on a surface element, itis computationally expensive to measure the illuminance at an arbitraryposition and orientation in the virtual environment, such as would beperformed with a photometer in a physical environment. The calculationsbecome prohibitively expensive when hundreds of meter positions must becalculated for climate-based annual daylight metrics such as spatialDaylight Autonomy (sDA) and Annual Sunlight Exposure (ASE), where themeters are typically positioned 0.6 meters above the floor.

In another novel contribution of the present invention, a virtualphotosensor is modeled as a transparent surface element with a finitearea. As previously described and in accordance with U.S. Pat. No.4,928,250, each iteration of the radiative flux transfer processtransferred unsent flux from a selected light source or surface elementto all other elements visible to that light source or element. However,each iteration then continues by transferring the same unsent flux toall photosensor elements visible to that light source or element. Thephotosensor elements record the flux they have received, but in beingtransparent they do not absorb or reflect any of the flux. Once theradiative flux transfer process has converged to a solution, thephotosensors have measured the incident spectral radiant flux and henceilluminance without participating in the radiative flux transferprocess.

Transition Surfaces

A disadvantage of representing diffuse daylight as n′ separateinfinite-distance light sources that emit parallel beams of light isthat the discrete nature of the beams becomes evident when projectedthrough small windows and openings separating the exterior environmentfrom the interior environment (FIG. 3). These parallel beams tend toproduce visible “spokes” of light patterns on the interior surfaces.Worse, they may contribute to errors in the subsequent radiative fluxtransfer calculations for the interior surface elements.

A prior art solution as implemented for example by AGi32 is to use avirtual camera to capture a hemispherical (“fisheye”) image of theexterior environment 300 as seen from the center of each window oropening 310 (referred to as a “transition” surface) or portion thereof(As will be understood by one skilled in the art, a virtual cameragenerates a two-dimensional image of the virtual three-dimensionalenvironment in the same manner as a digital camera would capture animage of the equivalent physical environment.) One method of captureutilizes a “hemicube” or similar geometric construct as disclosed inU.S. Pat. No. 4,928,250. This image is then projected from said centerinto the interior environment 320. Each pixel of the virtual imageprojects the average luminance or spectral radiant exitance of theexterior elements it represents onto the interior elements, therebyeffecting radiative flux transfer. (A portion of the luminous flux orspectral radiant flux may be absorbed or reflected by a window surface.)

It is evident that each pixel of the virtual image represents a ray ofluminous or spectral radiant flux emanating from the exteriorenvironment that intersects the window or opening at the position of thevirtual camera before transitioning into the interior environment. It isfurther evident that the angle of incidence θ between said ray and thewindow surface normal can be precalculated. The reduction in transmittedflux (that is, the specular transmission coefficient) through atransparent window (glass, plastic, or more generally a dielectricmaterial) due to Fresnel reflection can therefore be calculated on aper-pixel basis according to Schlick's approximation to the Fresnelequations:

T(θ)=1−R(θ)=1−(R ₀+(1−R ₀)(1−cos(θ))⁵)  (1)

where R₀ is the surface reflectance at normal incidence (i.e., θ=0).

In cases where the dimensions of the window or opening are large withrespect to the dimensions of the interior environment, it isadvantageous to spatially divide the window or opening into one or amultiplicity of subareas called “patches,” with the virtual camera(hemicube) positioned at the center of each patch.

A disadvantage of this approach is that the images project only luminousor spectral radiant flux. In order to project the luminous or spectralradiant flux distribution of the exterior environment into the interiorenvironment for each of the ‘n’ canonical radiosity solutions, it wouldbe necessary to repeat the transfer process ‘n’ times for eachtransition surface. This process involves calculating the form factorsof all of the exterior and interior environment elements visible fromthe transition surface, which as previously noted is computationallyvery expensive.

The present invention avoids this computational expense by firstassigning a unique identifier to every finite element in the exteriorand interior environments, then transferring the exterior environmentelement identifiers associated with each pixel rather than their averageluminance or spectral radiant exitance. When these pixels are“projected” onto the interior environment elements, it is possible toaccess the ‘n’ luminance or spectral radiant exitance values assigned toeach exterior environment element through their assigned identifiers. Inthis manner, the ‘n’ radiosity solutions for the exterior environmentcan be transferred to the interior environment without any significantcomputational burden.

The ‘n’ luminance or spectral radiant exitance values assigned to eachexterior environment element through their assigned identifiers may beconsidered attributes of the exterior environment elements. Without lossof generality, any set of element attributes may be accessed throughtheir assigned identifiers and transferring them to their correspondinginterior environment elements.

The process of accessing the ‘n’ luminance or spectral radiant exitancevalues assigned to each exterior environment element and transferringthem to its corresponding interior environment element iscomputationally expensive in that the volume of data typically exceedsthe cache size of current-generation CPUs. The inevitable cache missesresult in excessive memory bus traffic and poor CPU utilization,particularly for large environments with hundreds of transition surfacesand thousands of thousands of transition surface patches.

Image capture and projection utilizing a hemicube or similar geometricconstruct utilizes the well-known Z-buffering technique wherein thedistance from the virtual camera to each pixel (representing a point onthe surface of an opaque geometric element) is stored in atwo-dimensional “depth buffer” (also known as a “Z-buffer”), along withthe element identifier. Each opaque geometric element is scanned usingwell-known scanline techniques separately; if another surface pointscans to the same pixel, the distance of the closer element and itsidentifier are stored in the buffer pixel.

In a preferred embodiment, each transition surface is subdivided intoone or more patches, as disclosed herein. Assuming that a giventransition surface is comprised of ‘n’ patches, a depth buffer istemporarily assigned to the surface for the purpose of transferringelement identifiers, and its pixels randomly assigned integer “patchidentifiers” in the range of 1 to ‘n’.

Referring to FIG. 12, a transition surface comprised of two patches 1210and 1220 is shown wherein each patch is assigned a patch identifier.Virtual images 1230 and 1240 of an element 1200 in the exteriorenvironment are then captured with the virtual camera (not shown)positioned at transition surface patch centers 1215 and 1225respectively, As each exterior element is scanned, each of its pixels iscompared with the corresponding depth buffer pixel only if the patchidentifier matches that of the randomly assigned depth buffer pixel. Theexterior element 1200 is thereby represented in virtual images 1230 and1240 by the set of pixels 1235 and 1245 respectively.

After the set of ‘n’ transition surface patches have been processed, thedepth buffer will contain the equivalent of ‘n’ virtual images (one perpatch). Each pixel position is then assigned a random number between 1and ‘n’ such that each pixel position is unique to one of the ‘n’images. In FIG. 13, for example, two virtual images 1300 and 1310consisting of 64 pixels each have had their pixels randomly butexclusively chosen and combined in depth buffer image 1320.

These interspersed images are then projected onto the interiorenvironment from each of the ‘n’ patch centers. Referring, for example,to FIG. 14, an interspersed depth buffer image is projected from patchcenters 1405 and 1415 respectively of patches 1400 and 1410 ontoelements 1430, 1435, and 1440 of the interior environment.

Referring to FIG. 15, the method comprises: assigning a uniqueidentifier to each surface element in the virtual representation of theenvironment (Step 1510); subdividing each window or opening into one ora multiplicity of patches (Step 1512); capturing an image of visibleexterior surface elements as seen from the center of each patch of theone or more windows or openings, wherein the value of each pixel is theexterior element identifier visible to that pixel (Step 1514); assigningan empty image to each of the one or more windows or openings (Step1516); for each pixel of each empty image, selecting the correspondingpixel of a randomly selected image from the set of patch images for eachof the one or more windows or openings and assigning the exteriorelement identifier as the empty image pixel value (Step 1518); for eachof the one or more windows or openings, projecting its image onto thevisible surface elements of the interior elements of the environment asseen through the window or the opening (Step 1520); and indirectlydetermining luminance or spectral radiant exitance transferred from anexterior surface element to a corresponding interior surface element forone of a multiplicity of canonical radiosity solutions by way of theunique identifier of the exterior surface element (Step 1522).

In a preferred embodiment, the distance between the window patch centersand the nearest interior elements is greater than approximately twicethe distance between adjacent patch centers. If this condition cannot besatisfied, then it is advantageous to further subdivide the transitionsurface into smaller patches such that the condition is satisfied.)

As previously noted, each pixel of the virtual image represents a ray ofluminous or spectral radiant flux emanating from the exteriorenvironment that intersects the window or opening at the position of thevirtual camera before transitioning into the interior environment,irrespectively of whether the flux is represented directly or indirectlyby means of element identifiers. Suppose the image has ‘m’ pixels; bychoosing to randomly select pixels from ‘n’ virtual images for captureand projection, the above is equivalent to employing m/n rays for eachtransition surface patch.

The advantage of this embodiment is that the exterior elementidentifiers are projected through the transition surfaces onto theinterior environment once per transition surface rather than once pertransition surface patch. Compared to the process of accessing the ‘n’luminance or spectral radiant exitance values assigned to each exteriorenvironment element and transferring them to its corresponding interiorenvironment element, the process of capturing and projecting virtualimages is computationally inexpensive, especially if it is implementedin hardware using one or more graphics processor units (GPU). Thecomputational expense of accessing luminance or spectral radiantexitance values is thereby reduced by a factor of ‘p’, where ‘p’ is theaverage number of patches per transition surface for the environment.

This preferred embodiment is particularly useful for environments thatinclude a large number of glazed surfaces (or example, greenhouses). Itis also useful for calculating spatial Daylight Autonomy (sDA) andAnnual Sunlight Exposure (ASE) metrics, which require windows to besubdivided into patches no larger than 0.3 m×0.3 m.

A “half-space” may be defined as the three-dimensional space divided bythe plane of the transition surface that includes either the exteriorenvironment or the interior environment, while a “direction” is uniquelydefined by two orthogonal direction coordinates (for example, a verticalangle θ and a horizontal angle φ in a polar coordinates system) in thehalf-space, and a “position” is uniquely defined by two orthogonalposition coordinates (for example, ‘x’ and ‘y’ in a Cartesian coordinatesystem).

A virtual camera may be used that has a defined position on thetransition surface, while each ray corresponding to a pixel of a virtualimage has a defined direction within the exterior environment. Thus,each ray is uniquely defined by its coordinates in a four-dimensionalspace comprised of two direction coordinates and two positioncoordinates.

Each ray position may be chosen from anywhere on the transition surface.The position may be chosen from a grid of positions (corresponding, forexample, to the center of each transition surface patch), or randomlysuch that the set of positions is uniformly distributed across thetransition surface. (If the transition surface has been adaptivelysubdivided in accordance with the nearest exterior or interior elements,then the set of positions may be randomly chosen such that approximatelythe same number of positions are chosen from each patch.)

Each ray may be traced from its position on the transition surface inthe given direction to the closest exterior environment element,whereupon it is assigned the exterior environment identifier. The ray isthen traced from its position on the transition surface in the oppositedirection (that is, into the half-space defined by the transitionsurface and the interior environment) to the closest interiorenvironment element, whereupon the exterior environment elementidentifier as an attribute to the interior environment element.

Referring to FIG. 16, the method may comprise: assigning a uniqueidentifier to each element in the virtual representation of the exteriorand interior environments (Step 1610); selecting a set of directions inthe half-space defined by the transition surface and the exteriorenvironment (Step 1612); choosing a position on the transition surfacefor each direction (Step 1614); tracing a ray from the position on thetransition surface in the given direction to the closest exteriorenvironment element (Step 1616); assigning the exterior environmentelement identifier to the ray (Step 1618); tracing a ray from theposition on the transition surface in the opposite direction to theclosest interior environment element (Step 1620); assigning the exteriorenvironment element identifier as an attribute to the interiorenvironment element (Step 1622); and indirectly determining an elementattribute transferred from an exterior environment element to acorresponding interior environment element by way of the uniqueidentifier of the exterior environment element (Step 1624).

The above approaches assume that the windows are fully transparent.Recent advances in daylight control for high performance buildings havehowever introduced glazing systems with complex bidirectional scatteringdistribution functions (BSDFs), including light redirecting materialssuch as frits, ground glass, opal glass and plastics, prismatic glassand plastics, holographic films, diffraction gratings, and Fresnel lensarrays. The BSDFs may be spectrally neutral or (particularly withholographic and diffraction gratings) spectrally dependent.

In particular, the BSDFs of opal glass and plastics, which are comprisedof a dispersion of colloidal particles with a refractive indexsignificantly different from that of the matrix material, are determinedmostly by Tyndall scattering, wherein the degree of scattering is afunction of the fourth power of the wavelength of the incident light.Thus, while visible light with wavelengths of 400 nm to 700 nm may besubject to approximately the same degree of scattering, infraredradiation with wavelengths of 1 to 10 microns will be subject tosignificantly less scattering.

To address such systems, it is noted that such glazing systems bothredirect and scatter the incident flux. Thus, rather than simplyprojecting a captured virtual image into the interior environment, eachpixel can be assigned a scattering distribution function that describeshow an incident ray with given altitudinal and azimuthal anglescorresponding to the pixel will be scattered. These functions can becompactly represented using spherical and hemispherical harmonics,wavelets, radial basis functions, and other two-dimensional signalcompression techniques as will be known to those skilled in the art ofimage compression. The BSDF function coefficients may be measured orcalculated from virtual representations of the glazing systems.

Numerous bidirectional reflectance distribution function (BRDF) modelshave been developed for computer graphics applications. Their primaryintent is to generate photorealistic reflections from glossy surfaces incomputer games, architectural renderings, and computer graphicsanimations. For the purposes of daylight simulations, however, it isessential that the BRDF (or more generally, BSDF) model is energyconservative.

A BSDF model is energy conservative if no more light can be reflectedfrom a point on a surface than the amount of light incident on thatpoint. Many models include in their formulation a normalizationfactor—different for each model—that is necessary to ensure energyconservation as the model parameters are changed.

A widely-used BRDF model for representing glossy reflections from opaquematerials such as plastic is the Blinn-Phong model, and definedmathematically as:

L _(r)(i,r)=R(n·h)  (2)

where in FIG. 19, i is the incident light vector, r is the reflectedlight vector, R is the diffuse reflectance of surface, m is theBlinn-Phong exponent, n is the surface normal,

$h = \frac{i + r}{{i + r}}$

is the normalized half-vector between the incident and reflected lightvectors, and L_(r)(i,r) is the reflected luminance (assuming an incidentray with unit luminance). The reflection from the material surface isfully diffuse for m=0, but becomes increasingly specular with increasingm, and becomes fully specular for m=∞.

An equivalent formulation of the Blinn-Phong model is, again inreference to FIG. 19:

L _(r)(i,r)=R(cos θ_(n))^(m)  (3)

where θ_(n) is the angle between the surface normal n and the normalizedhalf-vector h.

Unfortunately, the Blinn-Phong BRDF model is not energy conservative. Inorder to ensure that the total amount of light remains constant as theBlinn-Phong exponent m is varied, it is necessary to add amultiplicative normalization term. A common formulation is:

$\begin{matrix}{{L_{r}\left( {i,r} \right)} = {\frac{m + 8}{8\pi}{R\left( {n \cdot h} \right)}^{m}}} & (4)\end{matrix}$

This is known, however, to be an approximation suitable for computergraphics applications only. As such, it is not suitable for use withdaylight simulation programs.

To determine the exact normalization term, the total amount of light Φreflected from the surface is determined by integrating the reflectedluminance L_(r)(i,r) over the hemisphere:

Φ=∫_(ω)(θ,φ)cos θdω  (5)

which, expressed in spherical coordinates, is:

Φ=∫₀ ^(2π)∫₀ ^(π/2) L(θ,φ)cos θ sin θdθdφ  (6)

If the incident light is perpendicular to the surface, the reflectedluminance L is independent of φ and we have θ being the angle betweenthe reflectance vector r and the surface normal n, which means thatθ_(n)=θ/2. Thus:

$\begin{matrix}{\Phi = {2\pi \; R{\int_{0}^{\pi/2}{\left( {\cos \left( \frac{\theta}{2} \right)} \right)^{m}\cos \; {\theta sin}\; \theta \; d\; \theta}}}} & (7)\end{matrix}$

Using the half-angle formula

${\cos \left( \frac{\theta}{2} \right)} = \sqrt{\left( {1 + {\cos \; \theta}} \right)/2}$

and the substitution t=cos θ (which gives dt=−sin θdθ), we have:

$\begin{matrix}{\Phi = {{2\; \pi \; {R\left( {- {\int_{1}^{0}\left( \sqrt{\frac{1 + t}{2}}\  \right)^{m}}} \right)}{tdt}} = {2\; \pi \; R{\int_{0}^{1}{\left( \frac{1 + t}{2}\  \right)^{m/2}{tdt}}}}}} & (8)\end{matrix}$

which can be evaluated using integration by parts (that is,

${\int{udv}} = {{uv} = {{- {\int{{vdu}\text{)}\mspace{14mu} {using}\mspace{14mu} u}}} = \left( \frac{1 + t}{2}\  \right)^{m/2}}}$

and v=t, and by removing the constant reflectance R, we have:

$\begin{matrix}{{2\; \pi \left( {\left\lbrack {\frac{4}{m + 2}{t\left( \frac{1 + t}{2}\  \right)}^{{({m + 2})}/2}} \right\rbrack_{t = 0}^{1} - {\frac{4}{m + 2}{\int_{0}^{1}{\left( \frac{1 + t}{2}\  \right)^{{({m + 2})}/2}{dt}}}}} \right)} = {{\frac{8\; \pi}{m + 2}\left( {\left\lbrack {t\left( \frac{1 + t}{2}\  \right)}^{{({m + 2})}/2} \right\rbrack_{t = 0}^{1} - {\frac{4}{m + 4}\left\lbrack \left( \frac{1 + t}{2}\  \right)^{{({m + 4})}/2} \right\rbrack}_{t = 0}^{1}} \right)} = {{\frac{8\; \pi}{m + 2}\left\lbrack {{\cos \; {{\theta cos}\left( \frac{\theta}{2} \right)}^{m + 2}} - {\frac{4}{m + 4}{\cos \left( \frac{\theta}{2} \right)}^{m + 4}}} \right\rbrack}_{\theta = {\pi/2}}^{0} = {\frac{8\; {\pi \left\lbrack {{\cos \left( \frac{\theta}{2} \right)}^{m + 2}\left( {{4\left( {\cos \left( \frac{\theta}{2} \right)} \right)^{2}} - {\left( {m + 4} \right)\cos \; \theta}} \right)} \right\rbrack}_{\theta = 0}^{\pi/2}}{\left( {m + 2} \right)\left( {m + 4} \right)} = {\frac{8\; {\pi \left\lbrack {\cos \; \left( \frac{\theta}{2} \right)^{m + 2}\left( {{2\left( {1 + {\cos \; \theta}} \right)} - {\left( {m + 4} \right)\cos \; \theta}} \right)} \right\rbrack}_{\theta = 0}^{\pi/2}}{\left( {m + 2} \right)\left( {m + 4} \right)} = {\frac{8\; {\pi \left\lbrack {\cos \; \left( \frac{\theta}{2} \right)^{m + 2}\left( {2 - {\left( {m + 2} \right)\cos \; \theta}} \right)} \right\rbrack}_{\theta = 0}^{\pi/2}}{\left( {m + 2} \right)\left( {m + 4} \right)} = {\frac{8\; {\pi \left\lbrack {2^{1 - {({{({m + 2})}/2})}} - \left( {2 - \left( {m + 2} \right)} \right)} \right\rbrack}}{\left( {m + 2} \right)\left( {m + 4} \right)} = \frac{8\; {\pi \left( {2^{{- m}/2} + m} \right)}}{\left( {m + 2} \right)\left( {m + 4} \right)}}}}}}}} & (9)\end{matrix}$

The exact energy-conserving normalization term is therefore

$\frac{\left( {m + 2} \right)\left( {m + 4} \right)}{8\; {\pi \left( {2^{{- m}/2} + m} \right)}},$

and so the exact energy-conserving Blinn-Phong BRDF model is:

$\begin{matrix}{{L_{r}\left( {i,r} \right)} = {\frac{\left( {m + 2} \right)\left( {m + 4} \right)}{8\; {\pi \left( {2^{{- m}/2} + m} \right)}}{R\left( {n \cdot h} \right)}^{m}}} & (10)\end{matrix}$

where both normalization terms are shown in FIG. 20.

Plotted against the approximate normalization term

$\frac{m + 8}{8\; \pi},$

the error is shown in FIG. 21. The maximum error for m=10 is 7.5percent, which justifies the use of the exact normalization term inphysically correct global illumination calculations.

Modeling BSDFs

Referring to FIG. 22, it is possible to model the BTDFs of translucentmaterials by adapting the Blinn-Phong BRDF model as follows:

L _(t)(i,t)=T(n·h′)^(m)  (11)

where t is the transmitted light vector, T is the surface diffusetransmittance,

$h^{\prime} = \frac{{- i} + t}{{{- i} + t}}$

is the normalized half-vector between the reversed incident ray i* andtransmitted light ray t, and L_(t)(i,t) is the transmitted luminance oflight ray 1230 (again assuming an incident ray with unit luminance).Thus:

$\begin{matrix}{{L_{t}\left( {i,t} \right)} = {\frac{\left( {m + 2} \right)\left( {m + 4} \right)}{8\; {\pi \left( {2^{{- m}/2} + m} \right)}}{T\left( {n \cdot h^{\prime}} \right)}^{m}}} & (12)\end{matrix}$

where the exact normalization term

$\frac{\left( {m + 2} \right)\left( {m + 4} \right)}{8\; {\pi \left( {2^{{- m}/2} + m} \right)}}$

can, of course, be precalculated.

Together, Equations 10 and 12 provide computationally efficient andphysically correct analytic equations to represent the BSDFs of manycommon fenestration devices, with example BTDFs shown in FIG. 23 forvarious incidence angles 8 and Blinn-Phong exponents m. These equationsare easily integrated into a ray tracing program, such as Radiance,where they may serve as probability distribution functions for thescattered light rays.

As may be readily understood, a database of fenestration devices andtheir associated Blinn-Phong exponents may be maintained in a computerprogram.

Radiosity and BSDFs

It is not immediately obvious how Equations 10 and 12 can be integratedinto a radiosity-based program for daylight modeling. After all, theclassic radiosity method is only capable of modeling diffuse reflectionand transmission of light from surfaces. To understand how BSDFs can bemodeled using radiosity, it is first necessary to review the methoditself

The radiosity method is based on the concept of transferring luminous orradiant flux (i.e., light) between finite area surfaces (called“patches”) rather than by randomly tracing light rays. Referring to FIG.24, at each step of the calculation process, flux Φ is diffusely emittedor reflected from source patch S and received by a receiver patch R thatis visible to the source patch. (In practice, there are typically manyreceiver patches, some of which may partially or fully obscure thevisibility of other receiver patches as seen from the source patch.)

Form Factors

The portion of the flux that is diffusely emitted by or reflected fromsource patch S and which is received by receiver patch R is determinedby the “form factor” F between the source and receiver patches. It is anextremely difficult and arduous mathematical problem to calculate theform factor between two patches of arbitrary shape and dimension, but ageometric approach known as Nusselt's analogy (FIG. 25) can be used todetermine the form factor F between an infinitesimally small sourcepatch dS and a finite area receiver patch R.

The outline of the receiver patch as seen from the source patch isprojected onto the surface of an imaginary unit radius hemisphere Hcentered on the source patch. This projection P is then projecteddownwards onto the base B of the hemisphere; the area A of the downwardsprojection divided by the area of the base (i.e., π) is the form factorF. That is:

F=A/π  (13)

Unfortunately, solving for the form factor F of a receiver patch ofarbitrary shape and dimension using Nusselt's analogy is equallydifficult, as it involves projective spherical geometry. However, giventwo infinitesimally small patches dS and dR separated by ray E withlength r as shown in FIG. 26, the differential solid angle dω subtendedby dR as seen from dS is:

dω=cos θ_(R) dA _(R) /r ²  (14)

where dA_(R) is the differential area of dR and θ_(R) is the anglebetween the ray E and the surface normal n_(R) of receiver element dR.If we project the intersection of this solid angle with the unithemisphere onto the hemisphere base, the area of the projection is:

dA=cos θ_(S) dw=cos θ_(S) cos θ_(R) dA _(R) /r ²  (15)

where θ_(S) is the angle between the ray and the surface normal n_(S) ofthe source patch.

Substituting Equation 15 into Equation 13 then gives:

$\begin{matrix}{{dF} = \frac{\cos \; \theta_{S}\cos \; \theta_{R}{dA}_{R}}{\pi \; r^{2}}} & (16)\end{matrix}$

for the differential form factor dF from dS to dR.

Hemicubes

If the hemisphere H of FIG. 25 is replaced with a hemicube C (i.e., halfof a cube) as shown in FIG. 27, the projection P of a receiver patch Ronto the hemicube involves planar rather than spherical surfaces, and sothe determination of the form factor F becomes tractable. If, forexample, each face of the hemicube is subdivided into a grid of squarecells, the form factor F becomes:

F≈ΣΔF _(covered)  (17)

where ΔF_(covered) represents the delta form factor of each hemicubeface cell that the receiver patch R projects onto. As the size of thecells decreases, the delta form factor converges to dF, and so theapproximation converges to the form factor F.

Referring to FIG. 28, a hemicube face cell C is oriented such that itsnormal n_(C) is at an angle θ_(C) to the ray E with length r between thehemicube center S and the center of the face cell. From Equation 15, thedelta form factor ΔF_(C) is therefore:

$\begin{matrix}{{\Delta \; F_{C}} = {\frac{\cos \; \theta_{S}\cos \; \theta_{C}}{\pi \; r^{2}}\Delta \; A_{C}}} & (18)\end{matrix}$

where ΔA_(C) is the area of the face cell (as a fraction of the fullface area), and the angle θ_(S) is determined by:

cos θ_(S) =E·n _(S) /r  (19)

where n_(S) is the hemicube plane's normal.

Cubic Tetrahedrons

A cubic tetrahedron (a three-sided pyramid with an apex angle of 90degrees for each side) is more computationally efficient for the purposeof calculating form factors. As shown in FIG. 29, the cubic tetrahedronis oriented in three-dimensional Cartesian coordinates (with axes u, v,and n) such that the vertices of its triangular base are positioned atv₀={1, 1, −2}, v₁={1, −2, 1}, and v₂={−2, 1, 1} and lie on the plane P.With this, the cubic tetrahedron base normal n_(S) is described by thevector

$\left\{ {\frac{1}{\sqrt{3}},\frac{1}{\sqrt{3}},\frac{1}{\sqrt{3}},} \right\}.$

Thus, referring to FIG. 28:

cos θ_(S)=(E _(u) +E _(v) +E _(n))/√{square root over (3)}r  (20)

For cells on the cubic tetrahedron face perpendicular to the v-axis, theangle θ_(C) is determined by:

cos θ_(C) =−E·n _(C) /r  (21)

where the face cell normal n_(C) is described by the vector {0, −1, 0}.Also, the face lies on the plane v=1. Therefore:

cos θ_(C) =E _(v) /r=1/r  (22)

The same result can be derived for the other two faces. Thus, for anycubic tetrahedron cell:

$\begin{matrix}{{\Delta \; F_{C}} = {\frac{E_{u} + E_{v} + E_{n}}{\pi \; r^{4}\sqrt{3}}\Delta \; A_{C}}} & (23)\end{matrix}$

However:

r ² =E _(u) ² +E _(v) ² +E _(n) ²  (24)

and for each face, one of E_(u), E_(v), or E_(n) will always be one.Therefore:

$\begin{matrix}{{\Delta \; F_{C}} = {\frac{u + v + 1}{\pi \; \left( {u^{2} + v^{2} + 1} \right)^{2}\sqrt{3}}\Delta \; A_{C}}} & (25)\end{matrix}$

where u and v range from 1 to −2.

These delta form factors can be precalculated and accessed through alookup table when receiver patches are projected onto the cubictetrahedron faces during radiosity calculations. It is important tonote, however, that form factors apply to the calculation of diffuseflux transfer only; they cannot be applied to non-diffuse flux transferbetween specular surfaces, or between glossy surfaces whose opticalproperties are described by bidirectional scattering distributionfunctions.

Modeling BSDFs

It is possible to consider Equations 10 and 12 as representing theintensity distribution of a point light source. This distribution willbe dependent on the incidence angle of the incident ray, but the analogystill holds. IES RP-16-10, Nomenclature and Definitions for IlluminatingEngineering [5], defines luminous and radiant intensity in a givendirection as:

I=dΦ/dω  (26)

where dΦ is the differential radiant flux contained in an elemental conewith differential solid angle dω that is oriented in the givendirection.

Referring to FIG. 30 then, the radiant flux from a point source S withan intensity distribution I(θ,φ) that is received by a cubic tetrahedronface cell C with surface area ΔA_(C) and is at a distance r and at anangle θ_(C) to the face cell normal n_(C) is:

Φ_(C) =f _(W) I(θ,φ)cos θ_(C) dω  (27)

where the differential solid angle dω is approximately:

dω≈cos θ_(C) ΔA _(C) /r ²  (28)

For each cell with variable coordinates u and v, and n=1, we have fromEquation 24:

r=√{square root over (u ² +v ²+1)}  (29)

and from Equation 22:

cos θ_(C)=1/r  (30)

Thus the radiant flux received by each cubic tetrahedron cell is:

ΔΦ_(C) ≈I(θ,φ)/A _(C)/(u ² +v ²+1)^(3/2)  (31)

and similarly, for coordinates pairs u, n and v, n. (A similarderivation is possible for hemicube face cells.) The termΔA_(C)/(u²+v²+1)^(3/2) is referred as the “diffusion parameter” of thecubic tetrahedron cell, or more generally of a ray originating from thecubic tetrahedron center in the direction of the center of the cell.Further, as previously noted, the term I(θ,φ) may be an analyticcalculation of the BSDF function such as the Blinn-Phong model, orderived from BSDF measurements of a physical fenestration device.

Thus, if the receiver patch R is projected on n cells of a cubictetrahedron, the flux Φ_(R) received by the patch from the point sourceis approximately:

Φ_(R)≈Σ_(n) I(θ,φ)ΔA _(C)/(u ² +v ²+1)^(3/2)  (32)

which again converges to Φ_(R) as the cubic tetrahedron cell sizedecreases.

Examples of direct sunlight and diffuse daylight transmitted through atranslucent surface onto an interior surface for different Blinn-Phongexponents m are shown in in FIGS. 31-33.

Referring to FIG. 34, the method may comprise: assigning transitiondiffusion characteristics (such as, for example, a Blinn-Phong exponent)based on the diffusion properties of a transition surface material (Step3410); assigning a unique identifier to each element in the virtualrepresentation of the exterior and interior environments (Step 3412);selecting a first set of directions in the half-space defined by thetransition surface and the exterior environment (Step 3414); selecting asecond set of directions in the half-space defined by the transitionsurface and the interior environment (Step 3416); choosing a position onthe transition surface for each exterior direction (Step 3418);selecting a first direction from the first set of directions (Step3420); tracing a ray from the position on the transition surface in thefirst direction to the closest exterior environment element (Step 3422);assigning the exterior environment element identifier to the ray (Step3424); selecting a second direction from the second set of directions(Step 3426); tracing a ray from the position on the transition surfacein the second direction to the closest interior environment element(Step 3428); assigning the exterior environment element identifier as anattribute to the interior environment element (Step 3430); andindirectly determining an element attribute transferred from an exteriorenvironment element to a corresponding interior environment element byway of the unique identifier of the exterior environment element (Step3432).

In Step 3412, the first set of directions may be defined, for example,by subdividing the sky dome into 256 triangular patches and choosing thebarycenter of each patch, or by calculating the annual solar path for agiven geographic location and choosing solar positions for selectedtimes and dates.

In Step 3414, the second set of directions may be defined, for example,by the centers of the face cells of a cubic tetrahedron or a hemicube.

In an embodiment, Steps 3420 through 3432 are iterated, eithersequentially or in parallel, for each direction of the first set ofdirections, while Steps 3420 through 3432 are similarly iterated, eithersequentially or in parallel, for each direction in the second set ofdirections.

Interior Environment

For interior environments, the light sources are the direct sunlight,diffuse daylight, and reflected light from exterior environment elementstransferred through the transition surfaces, plus the artificial lightemitted by the luminaires.

For daylight harvesting purposes, one or more luminaires are typicallyassigned to “zones” that can be independently dimmed or switched. Forthe purposes of illustration, there are assumed to be ‘p’ luminairezones in interior environment, where ‘p’ may be less than or equal tothe number of luminaires.

As with the exterior environment, the present invention assigns amultiplicity of ‘n’+‘p’ unsent flux parameters and ‘n’+‘p’ luminous orspectral radiant exitance parameters to each interior environmentelement, where ‘n’ is the number of divisions of the sky dome and solarpath. Radiative flux transfer techniques can again be employed tocalculate the distribution of luminous or spectral radiant flux due tointerreflections between surface elements in the interior environmentfor each sky and solar patch and each luminaire zone. (Similar to thesky dome and solar patches, each luminaire zone is assigned unitluminous or spectral radiant intensity.) The result is the generation of‘n’+‘p’ separate canonical radiosity solutions for the interiorenvironment. (The same approach can be applied, of course, to luminairezones illuminating the exterior environment.)

In one embodiment, the set of radiosity solutions comprises a list of‘n’ luminous exitance or spectral radiant exitance values for eachexterior geometric element and ‘n’+‘p’ luminous exitance or spectralradiant exitance values for each interior geometric element. In anotherembodiment, each such list is compressed using for example wavelet,radial basis function, or similar lossless or lossy signal compressiontechniques as will familiar to those skilled in the art of datacommunications. Further compression may be achieved through correlationsbetween the exitance lists of geometrically adjacent or close elements,as will be familiar to those skilled in the art of multichannel datacompression for electroencephalograms, electrocardiograms, and the like.

It is further noted that the number of finite elements ‘m’ representingan interior environment, each with an array ‘n’+‘p’ of unsent flux andluminous exitance values, represents an ‘m’×(‘n’+′p′) array. This arraycan be compressed using for example two-dimensional singular valuedecomposition (SVD) or eigenvector techniques such as are used for imagecompression, and an approximate representation of a portion of thematrix reconstructed as required with limited memory.

Spectral Radiosity

According to prior art radiosity techniques, the spectral powerdistributions of daylight and electric light sources are represented bya multiplicity of spectral bands. In a typical embodiment, threespectral bands are employed, representing red, green, and blue light,wherein the spectral bands are defined by CIE XYZ tristimulus valuesrather than wavelength ranges. A radiosity solution for the environmentis then calculated for each spectral band. (As will be understood bythose skilled in the art, there are an infinite number of band-limitedspectral power distributions that can have the same tristimulus values.)

Upon completion of the radiosity calculations, a virtual photometer maybe used to determine the approximate spectral radiant exitance at anyselected point on a surface element, as represented by the radiositysolution value for each spectral band, using prior art techniques.Similarly, the approximate spectral radiant exitance at an arbitraryposition and orientation in the virtual environment can be determined,again using prior art techniques.

A disadvantage of this technique is that the spectral powerdistributions of daylight and electric sources can only be approximatedby a finite number of spectral bands. With, for example, three spectralbands representing red, green, and blue light, the width of eachspectral band is approximately 100 nm. Unfortunately, some lightingsimulations require spectral resolutions of 10 nm or less. Inhorticultural lighting, for example, the growth of healthy lettuce isstrongly influenced by whether the plants are irradiated by green LEDswith peak wavelengths of 510 nm, 520 nm or 520 nm (Johkan, M. et al.2012. “Effect of Green Light Wavelength and Intensity onPhotomorphogenesis and Photosynthesis in Lactuca sativa,” Environmentaland Experimental Botany 75:128-133).

A prior art brute force technique is to employ as many spectral bands asnecessary to accurately represent the light source spectral powerdistribution and calculate a radiosity solution for each band. For 5 nmspectral resolution, for example, 60 spectral bands are required torepresent visible light with wavelengths from 400 nm to 700 nm.Unfortunately, this approach typically requires excessive amounts ofmemory.

A solution to this problem becomes evident when it is realized that thespectral transmittance distributions of transparent and translucentmaterials such as glass, and the spectral reflectance distributions ofopaque materials, rarely feature complex spectral features. An analysispresented by Fairman, H. S., and M. H. Brill, 2004, “The PrincipalComponents of Reflectances,” Color Research and Application29(2):104-110, for example, demonstrated that the measured spectralreflectance distributions of 3,534 object color samples of both naturaland synthetic materials are relatively smooth, and can be represented bythree to six principal components. FIG. 35 shows, for example, threeprincipal components V₁, V₂ and V₃, plus a mean component V₀; anyphysically plausible spectral power distribution can be represented as aweighted sum of these spectral power distributions.

With three principal components, it is possible to generate an excellentapproximation to a randomly chosen measured spectral power distributionusing six principal components, and a good approximation using onlythree principal components (FIGS. 36 and 37 respectively). As shown byFairman and Brill, it is possible to generate a physically plausiblespectral reflectance distribution for most natural and syntheticmaterials, knowing only their three CIE XYZ tristimulus values.

Given this, it is therefore possible to generate a physically plausiblespectral irradiance distribution for a virtual photometer measurementusing only three and six spectral bands. Referring to FIG. 38, if it isassumed that the spectral power distributions of the light sources inthe environment (both daylight and electric lighting) have equal energyspectral power distributions 3810, the virtual photometer measurement3850 represents the direct spectral radiant flux 3820 received from thelight sources and spectral radiant flux 3830 that is reflected and/ortransmitted interreflected by one or a multiplicity of surfaces in theenvironment. With each reflection or transmission, the surface color3840 modifies the spectral power distribution of the spectral radiantflux that contributes to the virtual photometer spectral irradiancemeasurement. Virtual photometer measurement 3850 is referred to as a“canonical” spectral irradiance, due to it being based on an illuminantwith an equal energy spectral power distribution.

Referring to FIG. 39, it therefore follows that if the light sourceshave a complex spectral power distribution 3930, the spectral irradiancedistribution measured by the virtual photometer is simply the canonicalspectral irradiance distribution measured by the virtual photometer withan equal energy spectral power distribution 3920 multiplied on aper-wavelength basis by the complex spectral power distribution 3910.

The flowchart shown in FIG. 40 therefore illustrates the process 4000for calculating the spectral irradiance received by a virtual photometerdue to direct and indirect spectral radiant flux from a light sourcewith a complex spectral power distribution. In Step 4010, a radiositysolution with three to six spectral bands is calculated for the lightsource, which may be daylight or electric lighting, in which the lightsource has an equal energy power distribution. In Step 4020, a virtualphotometer measures the canonical spectral irradiance due to spectralradiant flux received directly from the equal energy light source andspectral radiant flux that has been reflected and/or transmittedinterreflected by one or a multiplicity of surfaces in the environment.In Step 4030, the virtual photometer measurement, comprised of three CIEXYZ tristimulus values or, more generally, three to six principalcomponent values, is transformed in accordance to Fairman and Brill intoa spectral irradiance distribution. In Step 4040, the virtual photometercanonical spectral irradiance distribution is multiplied on aper-wavelength basis with the light source complex spectral powerdistribution to yield the spectral irradiance distribution that ismeasured by the virtual photometer.

If the environment includes a plurality of light sources with differentcomplex spectral power distributions, separate radiosity solutions mustbe calculated (Step 4010); these steps can be performed in parallel. Thefollowing Steps 4020, 4030 and 4040 are then performed for eachradiosity solution, following which the spectral irradiancedistributions are summed.

Optimal Luminous Environment

It is important to note that while the calculation of the ‘n’ canonicalradiosity solutions for the exterior environment and ‘n’+‘p’ canonicalradiosity solutions for the interior environment may require minutes tohours of computation time, these calculations need to be calculated onlyonce for a given environment. Thereafter, the radiosity solution for anysky luminance distribution and any dimming or switched state for theluminaire zones can be calculated as a weighted sum of the canonicalradiosity solutions, where the weights are the luminances of theassociated sky dome and solar patches and the luminous intensities ofthe luminaires for each zone. These calculations can be performed inmilliseconds.

A “sky condition” is uniquely determined according to the Perez (orsimilar) sky model by the time of day, the Julian date, the geographicalposition of the site in terms of longitude and latitude, and the directnormal irradiance and diffuse horizontal irradiance, which determinesthe spatial distribution of sky luminance and solar position. For thissky condition, there will be a range of luminaire zone dimmer or switchsettings and automated fenestration states that provide a comfortableluminous environment for the occupants of the interior environment andminimize energy consumption. FIG. 8 shows for example a building 800with a window 810 that admits daylight 830 from the Sun 820 and sky,wherein the amount of daylight is determined by the weather 840.Daylight sensor 850 senses the amount of daylight entering the buildingand artificial light 860 received from electric light source 870 andcommunicates the information to controller 880, which subsequentlyadjusts the amount of artificial light.

These are often competing goals. Occupants will be interested in havingrelatively constant illuminance of their workplaces, minimal glare fromthe windows, and infrequent changes in the automated fenestrationdevices (especially motorized blinds). Minimizing energy consumption mayhowever dictate frequent changes in workplace illuminance andfenestration devices, especially on partly cloudy days where theavailable daylight may change rapidly.

Assuming however that these competing goals can be satisfactorilyresolved, the building controller operation will be as shown in FIG. 4.In Step 400, the building controller 45 reads the input devices 10 to 40shown in FIG. 1. In Step 405, the building controller 45 calculates anoptimal luminous environment. This environment will be simulated bysuitable weightings of the precalculated canonical radiosity solutions.This step may include thermal and air flow calculations in response tosolar heat gain. In Step 410, the controller 45 sets the dimmers orswitches the circuit of the luminaire zones (i.e., the electric lightingdevices 55) and in accordance with the weights of the ‘p’ canonicalradiosity solutions for the luminaire zones. Similarly, the settings ofthe automated fenestration devices are set in accordance with theaverage transmittance of the windows (i.e., the transition surfaces),and the HVAC devices are set in accordance with the thermal and air flowcalculations.

Greenhouse lighting will not typically need to be dimmed during the day,as its primary purpose is to provide supplemental lighting necessary forplant growth and health. It is however necessary for the greenhousecontroller 145 to record the daylight measured by photosensors 100 atregular intervals and calculate the Daily Light Integral. If the valueof this metric is less than that required for the specific plant, thensupplemental lighting may be used to meet the plant's dailyrequirements. The greenhouse controller operation will otherwise be thesame as shown in FIG. 4, although the thermal and air flow calculationswill likely need to account for changes in humidity due to evaporation.

FIG. 5A shows the process of reading the input devices for buildingcontroller 45 in more detail, namely reading the daylight photosensors10 (Step 500), reading the occupancy or vacancy sensors 15 (Step 505),reading the timers 20 (Step 510), reading the personal lighting controls25 (Step 515), reading the utility smart power meters 30 (Step 520),reading the current weather report feed 35 (Step 525), and reading thetemperature sensors 40 (Step 530) if present. Of course, the sequenceshown in FIG. 5 is for illustration purposes only, and so is arbitrary.

FIG. 5B shows the process of reading the input devices in more detailfor greenhouse controller 145, namely reading the daylight photosensors100 (Step 535), reading the temperature sensors 105 (Step 540), readingthe timers 110 (Step 545), reading the humidity sensors 115 (Step 550),reading the wind speed and direction sensors 120 (Step 555), reading thesoil moisture sensors 125 (Step 560), reading the carbon dioxide and/oroxygen gas sensors 130 (Step 565), reading the utility smart power meter135 (Step 570), and reading the current weather report feed 140 (Step575). Of course, the sequence shown in FIG. 5B is for illustrationpurposes only, and so is arbitrary.

FIG. 6 shows the process of simulating the interior luminous environment(i.e., the distribution of luminous exitance over the interiorenvironment elements). In Step 600, the current time and date isdetermined. This uniquely determines the solar altitudinal and azimuthalangles. In Step 605, the current sky condition is determined from thedirect normal irradiance and diffuse horizontal irradiance, which may beobtained for example from the weather report feed 35 in FIG. 1A or 140in FIG. 1B. It may also however be inferred from the daylightphotosensor readings, as will be disclosed below.

The luminaire zone, automated fenestration, and environmental devicesettings are determined in Step 610. These settings will in general beconstrained by the input devices data, including the daylightphotosensors 10, occupancy sensors 15, timers 20, personal lightingcontrols 25, utility power meters 30, weather report feeds 35, andtemperature sensors 40 shown in FIG. 1A, and the daylight photosensors100, temperature sensors 105, timers 110, humidity sensors 115, windsensors 120, soil moisture sensors 125, gas sensors 130, utility powermeters 135, and weather report feeds 140 shown in FIG. 1B.

Further constraints may of course be imposed or implied by therequirements of external systems accessed through the communicationports 55 or 155, such as for example an HVAC system controller, anenergy storage system controller, a building automation system, a smartpower grid, and so forth.

The interior luminous environment and optionally thermal environment iscalculated in Step 620 using a weighted summation of the precalculatedcanonical radiosity solutions as previously disclosed.

The interior luminous and thermal environments are examined in Step 625to determine whether the goals of comfortable luminous and thermalenvironments for the occupants, or minimum requirements for plant growthand health, and the minimization of energy consumption have beenattained. If not, Steps 610 through 625 are repeated.

FIG. 7A shows the output process of the controller 45, which consists ofsetting the luminaire zone dimmers and switches in Step 700, setting theautomated fenestration device states in Step 705, and setting the HVACdevice states in Step 710. Not shown is communication with externalsystems with which the controller 45 may exchange information andcommands.

FIG. 7B shows the output process of the greenhouse controller 145, whichconsists of setting the supplemental lighting dimmers and switches inStep 715, setting the automated fenestration device states in Step 720,setting the humidity control device states in Step 725, setting the CO2distribution device states in Step 730, and setting the HVAC devicestates in Step 735. Not shown is communication with external systemswith which greenhouse controller 145 may exchange information andcommands.

Optimization

Referring to FIG. 1A, the controller 45 comprises a neuro-fuzzy orsimilar logic inference engine that can be trained using for exampledeep learning techniques to generate optimal luminaire zone dimmer andswitch settings, automated fenestration device states, and optimal HVACdevice settings for any given sky condition or time sequence thereof tosatisfy the goals of comfortable luminous and thermal environments forthe occupants and the minimization of annual energy consumption.

Referring to FIG. 1B, the greenhouse controller 145 similarly comprisesa neuro-fuzzy or similar logic inference engine that can be trainedusing for example deep learning techniques to generate optimal luminairezone dimmer and switch settings, automated fenestration device states,and optimal HVAC and other environmental settings for any given skycondition or time sequence thereof to satisfy the goals of minimalrequirements for plant growth and health, and the minimization of annualenergy consumption.

More generally, a neuro-fuzzy or similar logic inference engine can betrained to optimize for multiple goals, including but limited to powerconsumption, water consumption, fuel consumption, and variousenvironmental factors, plant requirements, animal requirements,aquaculture requirements, human requirements, and waste reduction,wherein the satisfaction of such goals may be subject to short-term andlong-term constraints, including but not limited to varying power, fueland water costs, plant fertilizer and animal feed costs, environmentalcompliance requirements, and waste handling and removal costs.

As an example of multiple goal satisfaction, a building or greenhousemay derive electrical power from a combination of utility power, windturbines, biomass generators, and electrical storage batteries, whereinsaid power may be required for supplemental electric lighting andheating. The most cost-effective combination will involve utility powercosts, biomass generator capacity, and battery discharge states, all ofwhich may vary continually. Training the controller logic inferenceengine therefore involves both historical and current weather data, aswell as data obtained from ongoing system operation.

A key requirement of any logic inference engine is an input data setfrom which the engine can be trained to generate suitable outputs. Giventhe ability to quickly simulate the interior luminous environment usingthe radiosity-based method disclosed herein, suitable training data setsmay include the Typical Meteorological Year (TMY) historical weatherdata files available for most geographic locations in North America (orequivalent weather data for other continents). These files providehour-by-hour direct solar irradiance and diffuse irradiance values thatcan be used directly with the present invention.

As will be appreciated by those skilled in the art, there are manydifferent artificial intelligence engines that can be trained togenerate suitable output values for a range of input values; theneuro-fuzzy logic engine is merely one embodiment.

With the ability to train the building controller 45 or greenhousecontroller 145 offline prior to deployment using a virtualrepresentation of the environment, it becomes possible to design adaylight harvesting system with predictable energy savings performance.This includes the ability to locate and orient the daylight photosensorswithin the virtual interior environment such that the neuro-fuzzy logicengine can learn to correlate their outputs with known virtual skyconditions. When the controller is used with physical daylightphotosensors within the physical interior environment, it is likely thatthe controller will be able to infer the correct sky condition strictlyfrom the photosensor outputs without the need for direct normal anddiffuse horizontal irradiance values from the weather report feed.

In one embodiment, the desired spatial illuminance distribution due tothe combination of daylight and electric lighting across surfaces andworkplanes in the virtual representation of the environment for eachtime and Julian date under various sky conditions is densely sampledwith virtual photosensors spaced for example 30 centimeters on center.The predicted virtual photosensor outputs are arranged in a rectangulararray. Using singular value decomposition (SVD) or eigenanalysis, thedominant singular vectors or eigenvectors of the array are determined.These vectors typically require much less storage space than the fullmatrix, and can later be used to reconstruct an approximation of thefull matrix using for example truncated SVD for determining correlationsbetween the virtual photosensor outputs and the spatial illuminancedistributions.

In another embodiment, the virtual photosensor outputs are correlateddirectly with the dominant singular vectors or eigenvectors of thearrays, thereby eliminating the need to reconstruct approximations ofthe full matrices.

With multiple daylight photosensors 10 or 100, possibly mounted in eachluminaire, it further becomes possible for the controller to performdata fusion by dynamically combining the photosensor outputs. That is,the controller may autonomously learn during its training with the TMYdata set that a given weighting of photosensor outputs reliably predictsthe correct sky condition during the morning, while another weighting isbetter during the afternoon. These weighted photosensor outputs cansimilarly be correlated with the corresponding spatial illuminancedistributions or their corresponding dominant singular vectors oreigenvectors during the training period.

With multiple gas concentration sensors 130, it further becomes possiblefor the controller to perform data fusion by dynamically combining thegas sensor outputs. That is, the controller may autonomously learnduring its training with the TMY data set that a given weighting of gassensor outputs reliably predicts the desired gas concentration for agiven sky condition during the morning, while another weighting isbetter during the afternoon. These weighted gas sensor outputs cansimilarly be correlated with the corresponding spatial illuminancedistributions or their corresponding dominant singular vectors oreigenvectors during the training period.

Once the building controller 45 or greenhouse controller 145 has beensuitably trained, the three-dimensional finite element representation ofthe exterior and interior environments may optionally be removed fromthe controller memory, as it may be sufficient during operation of thecontroller to determine the most likely spatial illuminance distributionaccording to its correlation with the physical photosensor outputs. Thiseliminates the need to recalculate the spatial illuminance distributionin real time according to the photosensor outputs.

If on the other hand the three-dimensional finite element representationof the exterior and interior environments is retained in the controllermemory or otherwise readily accessible to the controller, the radiativeflux transfer can be iteratively recalculated with changes to thematerial reflectances, window transmittances, lamp lumen outputs, andother model parameters at regular intervals (for example, once per week)to improve the correlations between the predicted and measuredphotosensor outputs, or to suggest more appropriate positions andorientations for the physical photosensors.

Similarly, the controller may autonomously learn during its trainingwith the TMY data set that the goals of a comfortable luminous andthermal environment for the occupants and the minimization of annualenergy consumption are best satisfied with a dynamic configuration ofluminaire zones. This may be particularly easy to accomplish if eachluminaire is assigned an Internet address or other unique identifiersuch that a single dimming command can be sent to any number ofluminaires.

Further, controller training need not be limited to the TMY orequivalent weather data set. As will be known to those skilled in theart, the controller can learn to adapt its behavior to changes in theenvironment, such as changing user preferences for illumination oroccupancy sensor events, or to changes in the requirements of externalsystems such as for example an HVAC controller or smart power grid.

In another embodiment, the three-dimensional finite elementrepresentation of the interior environment is augmented with virtualrepresentation of occupants (for example, office workers) whosestochastic behavior includes movement within the environment to simulateoccupancy sensor events, and whose user preferences dim or switch theluminaires at various locations (such as for example open officecubicles or workstations). These virtual occupants can be used tosimilarly train the controller offline prior to deployment, and todetermine an optimal placement of occupancy sensors.

Following deployment, the controller can similarly record occupancysensor events and personal lighting control operations to refine thestochastic behavior of the virtual occupants.

During the controller training period, problems such as lamp failures,lamp lumen depreciation, photosensor and occupancy sensor failures, andnetwork communication failures can be simulated, thereby enabling thecontroller to learn appropriate response strategies to minimize energyconsumption, optimize spatial illuminance distributions, and satisfyother predetermined goals.

It will further be evident to those skilled in the art that neuro-fuzzylogic and similar artificial intelligence engines can be trained torecognize temporal patterns. In particular, the controller can implementa Kalman filter or other construct to observe a minute-by-minutesequence of sky conditions and predict the sky condition and occupantbehavior (as determined for example by the occupancy sensor outputs andpersonal lighting control commands) for some period of time into thefuture. Based on this ongoing sequence of predictions, it can develop along-term strategy of daylight harvesting system settings that willsatisfy the goals of occupant comfort and annual energy savings.

The controller may also receive input data from and send output data toa plurality of buildings that may be geographically proximate or remote.The controller may thereby learn commonalities and differences inoperational requirements between a possibly disparate set of buildingsthat further refines its behavior to changes in the environment.

Finally, in yet another embodiment, the physical controller is simulatedin software for the purposes of training during the design phase. Theresultant data from the training period is then used to program thephysical controller prior to deployment.

Controller Topology

In one embodiment, the building controller 45 or greenhouse controller145 is implemented as a centralized system. In another embodiment, thecontroller is implemented as a decentralized and self-organizing systemwith computational intelligence embedded in each luminaire andassociated “Internet of Things” hardware devices connected via abidirectional communications network.

The luminaire network topology may be based on swarm intelligence, suchas for example particle swarm or ant colony optimization, where thegoals are for example desired spatial illuminance distributions, minimalenergy consumption, user preferences, minimal visual glare, and soforth.

Reinforcement learning algorithms can be employed to provideself-learning capabilities for the controller in achieving these goals.In one embodiment, the controller can simultaneously operate thephysical daylight harvesting system in real time and model theenvironment offline using “what-if” scenarios to improve controllerperformance over time, wherein these scenarios may be implemented usingevolutionary algorithms, including genetic algorithms.

FIG. 9 shows a plurality of sensors and data feeds 910 providingreal-time information to a controller 920. In Step 930, the controllerreads the input data, then in Step 940 performs a virtual simulation ofthe building environment based on the virtual simulation parameters 980.The results of this simulation are compared with the predicted sensordata in Step 950. If the difference between the predicted and measureddata exceeds predefined limits, the simulation parameters are updatedand Steps 940 and 950 are repeated; otherwise Step 970 outputs commandsto the external devices 990.

Aquaculture

In another embodiment, the environment is an aquaculture systemcomprised of tanks, ponds, or ocean and lake enclosures used to raisefish for commercial purposes, where overhead or underwater lighting maybe used to control feeding, reproduction, circadian rhythms, and otherspecies management factors (US 2014/0355254). Supplemental lighting maybe required, particularly during winter months. Given the often remotelocations of marine net pens, providing electrical power often entailsconsiderable costs.

In addition to the Perez sky model or similar mathematical model forpredicting sky conditions, a mathematical model representing thedistribution of light in turbid fluid media may be employed. Such modelsenable the prediction of illuminance and other photometric parameters atspecified depths and distances for specified or measured waterturbidity. Additional parameters such as tidal depth and fish orshellfish behavior models may also be employed.

Aspects of the present disclosure may be implemented as a methodperformed by a predictive daylight harvesting system. FIG. 10 is a flowdiagram depicting an example method 1000.

At 1010, the method includes obtaining one or more input valuesidentifying a target distribution of direct and interreflected radiationwithin an environment. Additionally, or alternatively at 1010, themethod includes obtaining one or more input values identifying ageographic location, a time of day, and a date. Input values may beobtained received from external sources, internally generated, orretrieved from memory. As an example, a target distribution of radiationmay be user-defined, whereas the time and date may refer to the currenttime and date, or a future time and date. The geographic location mayrefer to a geographic location of a building or other environment thatis to be controlled.

At 1012, the method includes calculating an expected sky condition forthe geographic location, the time of day, and the date based on amathematical model. In an example, the mathematical model incorporateshistorical weather data. Alternatively, or additionally, themathematical model was previously trained on or using historical weatherdata. Alternatively, or additionally, the method may include obtaining acurrent weather feed from a remote device over a communications network,and calculating the expected sky condition may be further based oninformation contained in the current weather feed.

At 1014, the method includes calculating a predicted distribution ofdirect and interreflected solar radiation within the environment basedon the expected sky condition. In at least some implementations,calculating the predicted distribution of direct and interreflectedsolar radiation within the environment based on the expected skycondition includes determining a transfer of solar radiation through oneor more windows or openings within a virtual representation of theenvironment by: (1) assigning a unique identifier to each surfaceelement in the virtual representation of the environment; (2) capturingan image of visible surface elements of an exterior of the environmentas seen through a window or an opening of the one or more windows oropenings; (3) projecting the image onto the visible surface elements ofan interior of the environment as seen through the window or theopening; and (4) indirectly determining luminance or spectral radiantexitance transferred from an exterior surface element to a correspondinginterior surface element for one of a multiplicity of canonicalradiosity solutions by way of the unique identifier of the exteriorsurface element.

At 1016, the method includes obtaining measurement data from one or morephotosensors measuring an initial distribution of direct andinterreflected radiation within the environment. The initialdistribution of direct and interreflected radiation may includeradiation from solar and/or electrical lighting sources. In at leastsome implementations, the measurement data obtained from the one or morephotosensors provides a less dense spatial sampling of the initialdistribution of direct and interreflected radiation within theenvironment than the predicted distribution of direct and interreflectedsolar radiation within the environment.

At 1018, the method includes determining, based on the measurement dataand the predicted distribution of direct and interreflected solarradiation, a target distribution of direct and interreflected artificialelectromagnetic radiation produced by electrical lighting to achieve thetarget distribution of direct and interreflected radiation within theenvironment. In at least some implementations, the target distributionof direct and interreflected artificial electromagnetic radiationproduced by the electrical lighting may be calculated by performingradiative flux transfer calculations with a virtual representation ofthe environment. The virtual representation may include geometry andmaterial properties of objects that influence the distribution ofluminous flux in the environment. As an example, a photosensor in thevirtual representation of the environment may comprise a transparentsurface element with a finite area that records radiant flux transferredfrom a selected light source or surface element at each iteration of theradiative flux transfer calculations.

At 1020, the method includes setting output parameters to one or moredevices controlling admittance of solar radiation into the environmentand/or controlling the electrical lighting to modify the initialdistribution to achieve the target distribution of direct andinterreflected radiation within the environment. The target distributionof direct and interreflected radiation within the environment mayinclude solar and/or artificial electromagnetic radiation components.

In at least some implementations, the method at 1020 may further includesetting the output parameters to maximize the solar radiation admittedinto the environment as a component of the target distribution of directand interreflected radiation while maintaining one or more additionalenvironmental parameters within a predefined range. As an example, theoutput parameters may be set to supplement the solar radiation admittedinto the environment with the target distribution of direct andinterreflected artificial electromagnetic radiation produced by theelectrical lighting. Non-limiting examples of environment parametersinclude temperature, humidity, light spectrum, light intensity, and/orother parameters described herein.

In at least some implementations, the method may further includecalculating a distribution of solar heat gain within the environmentbased on the measurement data and/or the predicted distribution ofdirect and interreflected solar radiation. One or more effects of thedistribution of solar heat gain on one or more environmental parametersmay be calculated. A target value for each environmental parameter maybe obtained. Setting the output parameters to one or more devices mayfurther varying a relative proportion of solar radiation to artificialelectromagnetic radiation within the target distribution of direct andinterreflected radiation to achieve the target value for each of the oneor more environmental parameters. Additional measurement data may bereceived from one or more additional sensors measuring the one or moreenvironmental parameters within the environment. Calculating the one ormore effects of the distribution of solar heat gain on the one or moreenvironmental parameters may be further based on the additionalmeasurement data received from the one or more additional sensors.

In at least some implementations, the method may further includereceiving one or more information feeds from external sources. Themeasurement data, the information feeds, and the output parameters maybe used to train an artificial intelligence engine to control the one ormore devices in response to the measurement data and information feeds.

In at least some implementations, the methods and processes, or portionsthereof described herein may be tied to a computing system that includesone or more computing devices. Such methods and processes may beimplemented as a computer-application program or service, anapplication-programming interface (API), a library, and/or othercomputer-program type.

FIG. 11 is a schematic diagram depicting an example computing system1100. Computing system 1100 is a non-limiting example of a controlsystem or controller that forms part of a predictive daylight harvestingsystem. Computing system 1100 can perform the methods and processes, orportions thereof described herein. Computing system 1100 is shown insimplified form. Computing system 1100 may take the form of one or morepersonal computers, server computers, mobile computing devices,electronic controller devices, and/or other computing devices.

Computing system 1100 includes a logic subsystem 1110 and a storagesubsystem 1112.

Computing system 1100 may further include an input subsystem 1114, anoutput subsystem 1116, a communication subsystem 1118, and/or othercomponents not shown in FIG. 11.

Logic subsystem 1110 may include one or more physical logic devicesconfigured to execute instructions. For example, the logic subsystem maybe configured to execute instructions that are part of one or moreapplications, services, programs, routines, libraries, objects,components, data structures, or other logical constructs. Suchinstructions may be implemented to perform a task, implement a datatype, transform the state of one or more components, achieve a technicaleffect, or otherwise arrive at a desired result.

The logic subsystem may include one or more processors (as an example ofphysical logic devices) configured to execute software instructions.Additionally, or alternatively, the logic subsystem may include one ormore hardware and/or firmware logic machines (as an example of physicallogic devices) configured to execute hardware or firmware instructions.Processors of the logic subsystem may be single-core or multi-core, andthe instructions executed thereon may be configured for sequential,parallel, and/or distributed processing. Individual components of thelogic subsystem may be distributed among two or more separate devices,which may be remotely located and/or configured for coordinatedprocessing. Aspects of the logic subsystem may be virtualized andexecuted by remotely accessible, networked computing devices configuredin a cloud-computing configuration.

Storage subsystem 1112 includes one or more physical, non-transitorymemory devices configured to hold instructions executable by the logicsubsystem in non-transitory form, to implement the methods and processesdescribed herein. When such methods and processes are implemented, thestate of storage subsystem 1112 may be transformed—e.g., to holddifferent data. Storage subsystem 1112 may include removable and/orbuilt-in devices. Storage subsystem 1112 may include optical memorydevices, semiconductor memory devices, and/or magnetic memory devices,among other suitable forms.

Storage subsystem 1112 may include volatile, nonvolatile, dynamic,static, read/write, read-only, random-access, sequential-access,location-addressable, file-addressable, and/or content-addressabledevices. Aspects of logic subsystem 1110 and storage subsystem 1112 maybe integrated together into one or more hardware-logic components. Whilestorage subsystem 1112 includes one or more physical devices, aspects ofthe instructions described herein alternatively may be propagated by acommunication medium (e.g., an electromagnetic signal, an opticalsignal, etc.) that is not necessarily held by a physical device for afinite duration.

The terms “module,” “program,” and “engine” may be used to describe anaspect of computing system 1100 implemented to perform a particularfunction. In some cases, a module, program, or engine may beinstantiated via logic subsystem 1110 executing instructions held bystorage subsystem 1112. It will be understood that different modules,programs, and/or engines may be instantiated from the same application,service, code block, object, library, routine, API, function, etc.Likewise, the same module, program, and/or engine may be instantiated bydifferent applications, services, code blocks, objects, routines, APIs,functions, etc. The terms “module,” “program,” and “engine” mayencompass individual or groups of executable files, data files,libraries, drivers, scripts, database records, etc. A “service”, as usedherein, may refer to an application program executable across multipleuser sessions. A service may be available to one or more systemcomponents, programs, and/or other services. In some implementations, aservice may run on one or more server-computing devices.

Input subsystem 1114 may include or interface with one or moreuser-input devices such as a keyboard, mouse, touch screen, microphoneetc. Input subsystem 1114 may include or interface with one or moresensor devices, such as previously described herein. Non-limitingexamples of such sensor devices include daylight illuminance sensors,photosynthetic photon flux density sensors, ultraviolet radiationsensors, infrared radiation sensors, spectroradiometric sensors,temperature sensors, humidity sensors, soil moisture sensors, carbondioxide sensors, oxygen sensors, wind speed/direction sensors, waterturbidity sensors, gas concentration sensors, and plant phytometricsensors. Input subsystem 1114 may additionally receive input data frominformation feeds such as weather report feeds, building automationsystems, energy storage systems, solar panels, electrical powerutilities, and smart power grids.

Output subsystem 1116 may include or interface with one or moreuser-output devices such as a graphical display device, touch screen,audio speakers, etc. Output subsystem 1116 may include or otherwiseinterface with one or more devices that control admittance of solarradiation into an environment and/or control electrical lighting, aswell as the various other control devices described herein. Otherexample devices that may form part of the output subsystem or mayotherwise interface with the output system include luminaire switchesand dimmers, moveable shades, electrochromic windows, and heating,ventilation, and air conditioning equipment. Such devices may bereferred to as external devices if they are not integrated with anelectronic controller.

Communication subsystem 1118 may be configured to communicatively couplecomputing system 1100 with one or more other devices. Communicationsubsystem 1100 may include wired and/or wireless communication devicescompatible with one or more different communication protocols. Asnon-limiting examples, the communication subsystem may be configured forcommunication via a wired or wireless wide area network, local areanetwork, and/or personal area network. In an example, the communicationsubsystem may allow computing system 1100 to send and/or receivemessages to and/or from other devices via a communications network.

As described herein, a variety of information in the form of data may bemeasured, collected, received, stored, retrieved from storage,processed, analyzed, organized, copied, reported, and/or transmitted inraw and/or processed forms. Data includes a set of one or more values(i.e., data values) of one or more parameters or variables. Such valuesmay be quantitate or qualitative in nature. Data may be represented byone or more physical quantities, attributes, or characteristics of oneor more signals or object states.

An object state refers to a physical state of a tangible, physicalobject, such as a device or machine. Within the context of a computingsystem or other electronic system, an object state may include a valueof a bit stored in a memory cell or other suitable bistable ormultistable electronic circuit (e.g., flip-flop or latch) of a memorydevice. As an example, a value of a bit may be defined by a high or lowphysical voltage value of a memory cell, corresponding to values of 1 or0 for the bit, respectively.

Data represented by one or more signals (i.e., data signals) may bepropagated by a communication medium, in the form of electrical signals,electromagnetic signals, optical signals, etc. Data signals may becommunicated over one or more wired and/or wireless communications linksor paths. Data signals may take the form of or form part of a modulatedsignal or a non-modulated signal. Data signals may be formatted orotherwise organized into one or more messages, streams, packets,datagrams, and/or frames as defined by one or more communicationsprotocols.

Data may be represented in a variety of digital and/or analog forms.Within the context of digital data, an object state or signal componentrepresenting an individual data unit may be observed or identified ashaving a discrete value of two or more discrete values. Within thecontext of analog data, an object state or signal component representingan individual data unit may be observed or identified as having a valuewithin a range of non-quantized values.

A collection of data may take the form of a set instructions that areexecutable by a machine to perform one or more operations. Suchinstructions may be referred to as machine-readable instructions thatdirect the machine to perform one or more operations. A set ofinstructions may take the form of software or a portion thereof (e.g., asoftware component). Software may include firmware, an operating system,an application program or other program type, a software plug-in, asoftware update, a software module, a software routine, or othersoftware component.

An organized collection of data may take the form of a database systemor other suitable data structure (e.g., an electronic file). A databasesystem includes one or more databases that define relationships andassociations between and among data objects. As an example, a dataobject (e.g., a user identifier) that includes a set of one or more datavalues may be associated with one or more other data objects (e.g., auser setting). A database system may be integrated with or form part ofa software component.

Data may include metadata that describes other data. Metadata describingthe structure of other data, such as a relationship or association ofdata objects in a database may be referred to as structural metadata.Metadata describing the content of other data may be referred to asguide metadata. A collection of data may include metadata and other datadescribed by that metadata.

The configurations and/or approaches described herein are exemplary innature, and specific implementations or examples are not to beconsidered in a limiting sense, because numerous variations arepossible. The specific methods or processes described herein mayrepresent one or more of any number of processing strategies. As such,various acts illustrated may be performed in the sequence illustrated,in other sequences, in parallel, or in some cases omitted. Likewise, theorder of the above-described processes may be changed. The subjectmatter of the present disclosure includes all novel and nonobviouscombinations and subcombinations of the various methods, processes,systems and configurations, and other features, functions, acts, and/orproperties disclosed herein, as well as any and all equivalents thereof.

The embodiments of the invention may be varied in many ways. Suchvariations are not to be regarded as a departure from the spirit andscope of the invention, and all such modifications as would be obviousto one skilled in the art are intended to be included within the scopeof the claims.

Dynamic Lighting Simulation

A multiplicity of n canonical radiosity solutions are calculated for agiven virtual environment, wherein each said radiosity solutionrepresents the sum of the direct and indirect illumination of allsurfaces due to light emitted by a single sky or solar patch, andwherein each environment geometric element 4100 is assigned an array4110 of n parameters representing the luminous or spectral radiantexitance of the element according to the n canonical radiosity solutions(FIG. 41).

Once the n canonical radiosity solutions have been calculated, theradiosity solution for any given sky luminance distribution for aspecified time and date can be calculated as a weighted sum of thecanonical radiosity solutions, where the weights are the luminances ofthe associated sky and solar patches as predicted by a sky luminancedistribution model.

Given that the sky and solar patches are light emitters, the concept maybe trivially extended to electric light sources, including physicalluminaires with measured luminous intensity distributions anddiffusely-emitting light sources such as organic light-emitting diode(OLED) panels, and wherein each canonical radiosity solution representsthe sum of the direct and indirect illumination of all surfaces due tolight emitted by all the light sources assigned to a lighting channel

Similarly, a specific radiosity solution can be calculated as a weightedsum of the canonical radiosity solutions, where the weights are theluminous flux outputs of the light sources as determined by the set oflighting channel intensity settings. For color-changing luminaires, thegamut of colors can be represented by red-green-blue (RGB) virtual lightsources whose intensities are independently controlled by three lightingchannels.

It is further evident that the set of n canonical radiosity solutionsmay include both those due to sky and solar (i.e., daylight) patches andelectric light sources, and specific radiosity solutions derivedaccordingly. Specific radiosity solutions may be calculated inmilliseconds, as opposed to minutes to hours for prior art approaches.

The same approach can also be applied to “virtual photometers,” whichare comprised of user-defined positions and orientations in the 3Denvironment where the incident illuminance is calculated for eachcanonical radiosity solution. Photometer readings for specific radiositysolutions can then be calculated according to the assigned weights ofthe lighting channels and daylight patches.

The real-time 3D interactive display of photometric analysis of specificradiosity solutions may include as architectural visualizations. Asshown by prior art implementations, however (e.g., the AGi32 lightingdesign and analysis software program by Lighting Analysts), graphicalapplication programming interfaces (APIs) such as OpenGL withappropriate hardware acceleration by commodity graphics processing units(GPUs) can provide real-time 3D interactive display capabilities—butonly for precalculated static radiosity solutions with predeterminedlight source intensities and sky and solar patch luminances. (Eachgeometric element is comprised of three or four vertices, whose RGBvalues are bilinearly interpolated from the element spectral radiantexitances as described in, for example, the book “Radiosity: AProgrammer's Perspective”. These vertex RGB values are then used forOpenGL rendering purposes.)

Texture Maps

The OpenGL API, which is continually evolving, introduced the OpenGLShading Language (GLSL) in 2006. This programming language providedsoftware developers with more direct control of the GPU graphicspipeline. One of the features of this language is the ability to assignand modify user-defined attributes to each vertex of the environmentsurfaces. (Microsoft similarly offers its proprietary High Level ShadingLanguage (HLSL) for Direct3D.)

This suggests an approach wherein a multiplicity of canonical radiositysolutions is stored on a per-vertex basis as vertex attributes. Eachattribute is a 32-bit value, which is sufficient to represent theappropriately-scaled spectral radiant exitance for a single canonicalradiosity solution. These attributes can then be weighted by theassociated lighting channel intensity or daylight patch luminance valuesin real time for display purposes.

The problem with this approach is that the minimum number of per-vertexattributes that an OpenGL implementation must support is 16. Further,some of these attributes may be allocated for vendor-specific functions,further limiting the number of available per-vertex attributes.

On the other hand, commodity graphic cards in desktop and laptopcomputers typically include 512 megabytes to several gigabytes ofrandom-access memory. This memory is most often used to store texturemaps, including color, stencil, and depth buffers, which are specific tothe rendered images.

In the present invention, this memory is repurposed to store the ncanonical radiosity solutions on a per-vertex basis. As shown in FIG.42, this comprises the steps of:

-   -   1. Step 4200—calculate n canonical radiosity solutions. Each        solution is comprised of a linear array of m=p+r floating point        values representing the spectral (red-green-blue) radiant        exitances of the p geometric elements and (optionally)        luminances or spectral radiant exitances of the r daylight        patches.    -   2. Step 4210—bilinearly interpolate the spectral radiant        exitance and (optionally) luminance values of the element and        patch vertices.    -   3. Step 4220—store the n linear arrays in sequence as a        rectangular texture map in the graphics card texture memory. The        maximum dimensions of a texture are vendor- and        product-specific, typically ranging from 1024 to 8192 pixels.        The element addresses of the n linear arrays therefore need to        be appropriately mapped to the row and column addresses of the        texture map. (If needed, the linear arrays can be stored across        multiple texture maps.)

A particular advantage of storing the canonical radiosity solutions asone or more texture maps is that the memory is accessible to the GPU viaits high-bandwidth memory bus, rather than needing to be marshalled bythe CPU from main memory and transferred through the much slowergraphics card bus.

Data Storage

Once the canonical radiosity solutions have been calculated, it isnecessary to store the data and the associated geometry and materialsinformation for later display. While there are many possible datastorage formats, one example suffices to demonstrate the principles.

For architectural visualization applications involving daylight, it isconvenient to separate the environment geometry into “interior” and“exterior” elements. One reason is that the average illuminance ofinterior (e.g., room) surfaces is typically orders of magnitude lessthan that of exterior surfaces. By applying different scaling factor tothe interior and exterior vertex colors, more aesthetically pleasingimages that include both types of surfaces can be produced.

In some lighting designs, the lighting channels may not be consecutivelynumbered. An “active channel map” is then required that maps eachuser-defined lighting channel number to the index of one of thecanonical radiosity solutions.

It is further convenient to represent objects in the virtual environmentas “instances,” wherein each instance is comprised of one or more“surfaces.”

The example binary file format therefore consists often sections:

1. Header

2. Opaque Exterior Surface Blocks List

3. Opaque Interior Surface Blocks List

4. Textured Exterior Surface Blocks List

5. Textured Interior Surface Blocks List

6. Placeholder Exterior Surface Blocks List

7. Placeholder Interior Surface Blocks List

8. Transparent Surface Blocks List

9. Virtual Photometers List

10. Texture Maps Header

11. Texture Maps List

where each “block” represents the triangular elements comprising aninstance surface.

Header

The file header consists of the following:

Offset Type Name Size Description Notes 0x0000 DWORD version 1 Version 10x0004 BYTE ident 16 UUID identifer 2 0x0014 float ec_scale 1 Exteriorcolor scaling factor 3 0x0018 float ic_scale 1 Interior color scalingfactor 3 0x001C float max_x 1 Maximum x-axis limit 0x0020 float max_y 1Maximum y-axis limit 0x0024 float max_z 1 Maximum z-axis limit 0x0028float min_x 1 Minimum x-axis limit 0x002C float min_y 1 Minimum y-axislimit 0x0030 float min_z 1 Minimum z-axis limit 0x0034 DWORD num_oesurf1 Number of opaque exterior surfaces 0x0038 DWORD num_oisurf 1 Number ofopaque interior surfaces 0x003C DWORD num_oeind 1 Number of opaqueexterior surface indices 0x0040 DWORD num_oiind 1 Number of opaqueinterior surface indices 0x0044 DWORD num_oevert 1 Number of opaqueexterior surface vertices 0x0048 DWORD num_oivert 1 Number of opaqueinterior surface vertices 0x004C DWORD num_tesurf 1 Number of texturedexterior surfaces 0x0050 DWORD num_tisurf 1 Number of textured interiorsurfaces 0x0054 DWORD num_teind 1 Number of textured exterior surfaceindices 0x0058 DWORD num_tiind 1 Number of textured interior surfaceindices 0x005C DWORD num_tevert 1 Number of textured exterior surfacevertices 0x0060 DWORD num_tivert 1 Number of textured interior surfacevertices 0x0064 DWORD num_pesurf 1 Number of placeholder exteriorsurfaces 0x0068 DWORD num_pisurf 1 Number of placeholder interiorsurfaces 0x006C DWORD num_peind 1 Number of placeholder exterior surfaceindices 0x0070 DWORD num_piind 1 Number of placeholder interior surfaceindices 0x0074 DWORD num_pevert 1 Number of placeholder exterior surfacevertices 0x0078 DWORD num_pivert 1 Number of placeholder interiorsurface vertices 0x007C DWORD num_dsurf 1 Number of transparent surfaces0x0080 DWORD num_dindex 1 Number of transparent surface indices 0x0084DWORD num_dvertex 1 Number of transparent surface vertices 0x0088 longoeblist_offset 1 Opaque exterior surface blocks list offset 0x008C longoiblist_offset 1 Opaque interior surface blocks list offset 0x0090 longteblist_offset 1 Textured exterior surface blocks list offset 0x0094long tiblist_offset 1 Textured interior surface blocks list offset0x0098 long peblist_offset 1 Placeholder exterior surface blocks listoffset 0x009C long piblist_offset 1 Placeholder interior surface blockslist offset 0x00A0 long dblist_offset 1 Transparent surface blocks listoffset 0x00A4 long pmlist_offset 1 Photometer list offset 0x00A8 longtmlist_offset 1 Texture maps list offset 0x00AC WORD num_alc 1 Number ofactive lighting channels 4 0x00AE WORD num_meters 1 Number of virtualphotometers 0x00B0 WORD num_tex 1 Number of textures 0x00B2 WORD padding1 Reserved 0x00B4 short alc_ids 65 Active lighting channel identifiersarray 5, 6 0x0136 short chan_map 65 Channel map array 7, 8 0x01B8 floatcie_cat02_00 1 CIE CAT02[0][0] 0x01BC float cie_cat02_01 1 CIECAT02[0][1] 0x01C0 float cie_cat02_02 1 CIE CAT02[0][2] 0x01C4 floatcie_cat02_10 1 CIE CAT02[1][0] 0x01C8 float cie_cat02_11 1 CIECAT02[1][1] 0x01CC float cie_cat02_12 1 CIE CAT02[1][2] 0x01D0 floatcie_cat02_20 1 CIE CAT02[2][0] 0x01D4 float cie_cat02_21 1 CIECAT02[2][1] 0x01D8 float cie_cat02_22 1 CIE CAT02[2][1] 0x01DC BOOLday_flag 1 Daylight flag  9, 10

The daylight file header data consists of the following:

Offset Type Name Size Description Notes 0x01E0 DWORD units 1 Measurementunits 10 0x01E4 float compass 1 Compass rotation (degrees) 0x01E8 floatlatitude 1 Site latitude (degrees) 0x01EC float longitude 1 Sitelongitude (degrees) 0x01F0 WORD num_meters 1 Number of meters 0x01F2BYTE reserved 526 Reserved

Notes

1. Vertex RGB color values should be multiplied by ec_scale fordisplaying exterior surfaces and ic_scale for displaying interiorsurfaces.

2. The UUID (Universally Unique Identifier) provides a unique identifierfor the file.

3. Vertex RGB color values are not gamma-corrected.

4. For the purposes of num_alc, color-changing luminaircs are consideredto have a single active lighting channel with three independentsubchannels.

5. Each active lighting channel identifiers array (alc_ids) element canhave one of five values:

-   -   −l —Inactive channel    -   0—Red color channel only    -   1—Green color channel only    -   2—Blue color channel only    -   3—All color channels

6. The zeroth element of alc_ids represents constant intensityluminaires, and has a value of 3 (all color channels).

7. The channel map array (chan_map) elements map active lighting channelidentifiers to the vertex color offsets in the block color_table arrays(see below ).

8. The zeroth element of chan_map represents constant intensityluminaires, and has a value of 0.

9. The optional daylight parameters are required only if day_flag, isTRUE.

10. Valid units values are:

-   -   0—feet    -   1—meters    -   2—millimeters

Surface Blocks

OpenGL indexed vertex arrays require the same OpenGL state for eacharray. There are therefore four surface block types:

1. Opaque surfaces (exterior and interior);

2. Textured surfaces (exterior and interior);

3. Placeholder surfaces (exterior and interior); and

4. Transparent surfaces.

Opaque surfaces do not have assigned textures.

Textured surfaces are opaque surfaces that have assigned BGR texturemaps.

Placeholder surfaces are semitransparent surfaces that have assignedBGRA texture maps.

Transparent surfaces comprise interior and exterior glass and window(exterior-interior transition) surfaces.

Each surface block consists of the following:

Offset Type Name Size Description Notes 0x00 char entity_name 33  Entityname 1 0x21 BYTE layer_id 1 Surface layer identifier 3 0x22 WORD tex_id1 Texture identifier 4 0x24 WORD tex_param 1 Texture parameters 0x26WORD — 1 Padding 0x28 float avg_red 1 Average texture map redreflectance 0x2C float avg_green 1 Average texture map green reflectance0x30 float avg_blue 1 Average texture map blue reflectance 0x34 DWORDnum_index 1 Number of indices 0x38 DWORD num_vert 1 Number of vertices0x3C DWORD index_table num_index Vertex indices table 0x40 floatposition_table num_vert * 3 Vertex position table — float chan_tablenum_vert * q Vertex channel table 5, 6 — float day_table num_vert * dVertex daylight patch table 5, 7 — float coord_table num_vert * 2Texture coords table 8, 9 Notes 1. The null-terminated entity_name has amaximum length of 32 characters. This field can be used by theapplication program to identify and selectively display entities. 2.This field can be used by the application program to selectively displaysurface layers. 3. Texture identifiers are one-based; a tex_id of zeroindicates that no texture is associated with the surface. 4. For opaqueand textured surface (including placeholder) blocks, q = (num_alc + 1);for transparent surface blocks, q = 1. (The vertex colors of transparentand placeholder surfaces are dependent only on the material color.) 5.Vertex channel and daylight patch spectral radiant exitances are writtenin RGB order using the Ward LogLuv encoding format (Larson, G. W. 1998.“The LogLuv Encoding for Full Gamut, High Dynamic Ranges Images,”Journal of Graphics Tools 3(1): 15-31) with adaptive run-length encoding(Glassner, A. S. 1991. “Adaptive Run-Length Encoding,” in Graphics GemsII, pp. 89-92. San Diego, CA: Academic Press). 6. If num_alc > 0, opaqueand textured surface (including placeholder) vertex spectral radiantexitances for successive active lighting channels are stored in thechannel_table. For example, if num_alc = 2 and the active lightingchannels are 5 (color-changing luminaire) and 8 (white light luminaire),the first vertex would be represented as: Channel Red Green Blue Nonecolor_table[0] color_table[1] color_table[2] 5 color_table[3]color_table[4] color_table[5] 8 color_table6] color_table[7]color_table[8] 7. For opaque and textured surface blocks, d = 257 ifday_flag in the file header is TRUE; otherwise d = 0. For transparentsurface blocks, d = 0. 8. Texture parameters are as required by OpenGL.9. The vertex texture coordinates table is present only if tex_id isnon-zero.

Virtual Photometers List

Each of the num_meter records in the virtual photometers list iscomprised of:

Offset Type Name Size Description Notes 0x00 float xpos 1 Positionx-axis coordinate 0x04 float ypos 1 Position y-axis coordinate 0x08float zpos 1 Position z-axis coordinate 0x0C float xori 1 Orient x-axiscoordinate 0x10 float pori 1 Orient y-axis coordinate 0x14 float zori 1Orient z-axis coordinate 0x18 float chan_table q Channel values table 1— float day_table d Daylight patch values table 2 Notes 1. q =num_alc + 1. 2. d = 257 if day_flag in the file header is TRUE;otherwise d = 0.

Texture Maps

BGR (24-bit) and BGRA (32-bit) texture maps are stored in orderaccording to their (one-based) texture identifiers.

Each texture map is preceded by its texture map header:

Offset Type Name Size Description Notes 0x00 BOOL alpha_flag 1 Alphaflag 1 0x01 WORD num_row 1 Number of rows 0x03 WORD num_col 1 Number ofcolumns 0x05 long tm_offset 1 Texture map offset Notes 1. The texturemap type is BGRA if alpha_flag is TRUE; otherwise it is BGR. 2. Texturemaps do not have gamma correction applied.

OpenGL Rendering

To generate real-time dynamic lighting simulations using OpenGL orDirect3D, the geometric information is read from the data storage mediumand used to generate index vertex arrays, while the texture maps areread into the GPU memory. Similarly, the per-vertex channel and daylightpatch spectral radiant exitances are decompressed and used to generateone or more texture maps to represent the canonical radiosity solutions.

To display the real-time images (typically at 30 or more frames persecond), the CPU portion of the program manages the display and userinteraction of lighting channel settings and (optionally) the daylightpatch luminances according to the user-defined sky model. This minimalamount of information is transferred to the GPU for each frame, where itis stored in a “uniform” (i.e., global) data array. (In an alternativeembodiment, the setting and luminance values are transferred only ifthey have changed since the previous frame.)

In a first embodiment, the lighting channel settings are obtained fromuser controls, such as slider bars, that are displayed by the program.The user controls may optionally implement industry-standard dimmingcurves such as, for example, square law or sigmoid (“S”) curves.

In a second embodiment, the lighting channel settings are obtained bymeans of executing a predefined lighting script, whose instructionsdetermine the channel on/off settings and dimmer fade rates according toa time schedule. Such a script may, of course, include audio tracksynchronization cues and virtual camera path instructions.

In a third embodiment (shown in FIG. 43), the lighting channel settingsare derived from an external light console 4300 that generates a datastream 4310 of dimmer instructions based on, for example, the USITTDMX512-A communication protocol. These instructions directly control theintensities of the light sources displayed in the virtual environmentgenerated on the monitor of computer 4320.

GPU Graphics Pipeline

The GPU portion of the program, as implemented using the GLSL or HLSLshading languages, comprises two separate subprograms referred to as the“vertex shader” and the “fragment shader,” representing two stages ofthe GPU's graphics pipeline.

The vertex shader program is executed on a per-vertex basis for theentire environment, where the GPU has dozens to thousands of processorcores performing the program calculations in parallel. Expressed inpseudocode, this program is:

Get vertex constant color Initialize channel color texel array index IFopaque or textured surface FOR each active color channel IF channel isenabled Get channel color array texel Multiply channel color by channeldimmer setting Add channel color to vertex color Increment channel colortexel array index ENDIF ENDFOR ENDIF Perform inverse gamma correctionSet vertex texture coordinates Calculate vertex positionwith the equivalent flowchart shown in FIG. 44.

In Step 4400, the program fetches the vertex constant color from thetexture map and uses its value to initialize the display vertex color.Of the n canonical radiosity solutions, one solution may represent thedirect and indirect illuminance distribution within the environment forconstant (i.e., non-dimmable) electric light sources. Each vertextherefore has a constant (albeit possibly zero) color.

In Step 4405, the program initializes the vertex channel color texelarray index in order to access the spectral radiant exitances arrayassociated with the vertex. (A “texel” is the OpenGL terminology for atexture map pixel.)

In Step 4410, the program determines whether the vertex belongs to anopaque or a textured surface. If true, the program proceeds to Step4415; otherwise, it proceeds to Step 4445.

In Step 4415, the program iteratively determines whether there are moreactive color channels to be processed. If true, the program proceeds toStep 4420; otherwise, it proceeds to to Step 4445.

In Step 4420, the program determines whether the color channel iscurrently enabled. If, for example, the lighting channel associated withthe color channel has been disabled, the color channel will also bedisabled. If true, the program proceeds to Step 4425; otherwise, itproceeds to to Step 4445.

In Step 4425, the program fetches the vertex channel color texel fromthe texture map.

In Step 4430, the program multiplies the vertex channel colorrepresented by the texel with the associated global (or, in OpenGLterminology, “uniform”) channel dimmer setting or sky patch luminance.

In Step 4435, the vertex channel color is added to the display vertexcolor.

In Step 4440, the vertex channel color texel array index is incremented,following which program control returns to Step 4415 to process the nextactive color channel (if any).

In Step 4445, inverse gamma correction is applied to the display vertexcolor. For historical reasons related to cathode ray tube (CRT)technology, most computer displays exhibit a nonlinear relationshipbetween specified pixel intensity I (typically in the range of 0 to 255)and pixel luminance L, described by the relationship L=I²², where thenominal exponent 2.2 is referred to as the display gamma. To ensure thatthere is a linear relationship between the vertex color and thedisplayed pixel luminance, each component C of the red-green-bluetriplet must undergo an inverse gamma correction: C′=C^(0.45).

In Step 4450, the vertex texture coordinates are set as required for anyvertex shader program.

In Step 4455, the vertex position is calculated as required for anyvertex shader program.

The fragment shader program is executed on a per-pixel basis for theentire displayed image, where again the GPU has dozens to thousands ofprocessor cores performing the program calculations in parallel.Expressed in pseudocode, this program is:

IF textures enabled frag_color = color * texel_color ELSE frag_color =color ENDIF IF NOT transparent surface frag_color = frag_color *exposure IF grayscale display gsv = color.red * 0.2125 + color.green *0.7154 + color.blue * 0.0721 frag_color = gsv ELSE IF pseudocolordisplay gsv = color.red * 0.2125 + color.green * 0.7154 + color.blue *0.0721 IF gsv < 0.25 frag_color.red = 0.0 frag_color.green = 4.0 * gsvfrag_color.blue = 1.0 ELSE IF gsv < 0.5 frag_color.red = 0.0frag_color.green = 1.0 frag_color.blue = 2.0 − 4.0 * gsv ELSE IF gsv <0.75 frag_color.red = 4.0 * gsv − 2.0 frag_color.green = 1.0frag_color.blue = 0.0 ELSE frag_color.red = 1.0 frag_color.green = 4.0 −4.0 * gsv frag_color.blue = 1.0 ENDIF ENDIF ENDIFwith the equivalent flowchart shown in FIG. 5.

In Step 4500, the program determines whether textures are enabled forthe displayed image. If true, the program proceeds to Step 4505;otherwise, it proceeds to Step 4510.

In Step 4505, the fragment color is initialized as the input color timesthe texel color, following which program control proceeds to Step 4515.(The terms “fragment,” “input,” and “texel” are understood as being GLSLshading language terminology for fragment shaders.)

In Step 4510, the fragment color is initialized as the input color,following which program control proceeds to Step 4515.

In Step 4515, the program determines whether the surface is atransparent surface. If false, the program proceeds to Step 4545;otherwise, it proceeds to Step 4520.

The user may specify an exposure setting to lighten or darken thedisplayed image. In Step 4520, the fragment color is multiplied by theexposure setting (which by default is unity).

The user-specified display type may be RGB color, grayscale, orpseudocolor. In Step 4525, the program determines whether a grayscaledisplay has been specified. If true, the program proceeds to Step 4530;otherwise, it proceeds to Step 4535.

In Step 4530, the program calculates the fragment color as the luminanceL (i.e., its grayscale value) of its red-green-blue (RGB) values asL=0.2125*R+0.7154*G+0.0721*B, following which program control proceedsto Step 4545. (The constants assume a video display calibrated inaccordance with the ITU-R BT.709 standard for high-definition televisiondisplays; display calibrated in accordance with other standards willrequire different constants.)

In Step 4535, the program determines whether a pseudocolor display hasbeen specified. If true, the program proceeds to Step 4540; otherwise,it proceeds to Step 4545.

In Step 4540, the program calculates the fragment color as apseudocolor. As an example, an algorithm to calculate the pseudocolor ofa luminance value L in the range of 0.0 to 1.0 is:

IF L < 0.25 R = 0.0 G = 4.0 * L B = 1.0 ELSE IF L < 0.5 R = 0.0 G = 1.0B = 2.0 − 4.0 * L ELSE IF L < 0.75 R = 4.0 * L − 2.0 G = 1.0 B = 0.0ELSE R = 1.0 G = 4.0 − 4.0 * L B = 0.0 ENDIF

In Step 4545, the program output the fragment color value as requiredfor any fragment shader program.

Implemented as a combination of CPU and GPU procedures in a softwareprogram as disclosed herein, it is possible to interactively displaydynamic lighting of virtual environments with hundreds of thousands ofvertices and hundreds of lighting channels in real time on commoditydesktop and laptop computers.

This invention has been disclosed in terms specific to the OpenGLgraphical application programming interface, but it is equallyapplicable Direct3D and other graphical APIs that provide direct controlof the GPU graphics pipeline.

Dynamic Lighting Simulation Optimization

A neuro-fuzzy or similar logic inference engine that can be trainedusing for example deep learning techniques to generate optimal memoryand data storage including texture maps, data storage, headers, surfaceblocks, virtual photometers, rendering, CPU usage, and GPU usage fordynamic lighting simulation.

We claim:
 1. A method performed by a lighting modelling system fordisplaying real-time dynamic lighting simulations, the methodcomprising: encoding a texture map as a multiplicity of canonicalradiosity solutions, each representing a lighting channel; storing thesolutions in the texture memory of a graphics processing unit; generatea multiplicity of lighting channel intensity settings; store thelighting channel intensity settings; access the canonical radiositysolutions on a per-vertex basis with a vertex shader; multiply theassociated vertex channel colors by the lighting channel intensitysettings; and sum the resultant colors to generate a vertex color fordisplay.
 2. A method of claim 1 wherein the values of the vertexconstant color from the texture map are used to initialize the displayvertex color with one of the n canonical radiosity solutions used torepresent the direct and indirect illuminance distribution within theenvironment for constant electric light sources for each vertex.
 3. Amethod of claim 2 wherein the vertex channel color texel array index isinitialized in order to access the spectral radiant exitances arrayassociated with the vertex.
 4. A method of claim 3, wherein for one ormore active color channels to be processed: the vertex channel colortexel from the texture map is fetched; the vertex channel colorrepresented by the texel is multiplied with the associated globalchannel dimmer setting or sky patch luminance; the vertex channel coloris added to the display vertex color. the vertex channel color texelarray index is incremented, the next active color channel (if any) isprocessed.
 5. A method of claim 2 wherein the inverse gamma correctionis applied to the display vertex color, comprising: each component ofthe red-green-blue triplet undergoes an inverse gamma correction; thevertex texture coordinates are set as required for any vertex shaderprogram; the vertex position is calculated as required for any vertexshader program; and The fragment shader program is executed on aper-pixel basis for the entire displayed image, utilizing a multiplicityof GPU processor cores.
 6. A method of claim 3 wherein the inverse gammacorrection is applied to the display vertex color comprising: eachcomponent of the red-green-blue triplet undergoes an inverse gammacorrection; the vertex texture coordinates are set as required for anyvertex shader program; the vertex position is calculated as required forany vertex shader program; and The fragment shader program is executedon a per-pixel basis for the entire displayed image, utilizing amultiplicity of GPU processor cores.
 7. A method of claim 1, whereintextures are enabled for the displayed image and the fragment color isinitialized as the input color times the texel color.
 8. A method ofclaim 1, wherein the fragment color is initialized as the input color.9. A method of claim 7, wherein the surface is an opaque surfacecomprising the steps of: specifying an exposure setting to lighten ordarken the displayed image; the fragment color is multiplied by theexposure setting;
 10. A method of claim 8, wherein the surface is atextured surface comprising the steps of: specifying an exposure settingto lighten or darken the displayed image; the fragment color ismultiplied by the exposure setting;
 11. A method of claim 9, wherein agrayscale display has been specified comprising the steps of:calculating the fragment color as the luminance of its red-green-bluevalues; calculating the fragment color as a pseudocolor; and output thefragment color value as required for any fragment shader program.
 12. Amethod of claim 10, wherein a grayscale display has been specifiedcomprising the steps of: calculating the fragment color as the luminanceof its red-green-blue values; calculating the fragment color as apseudocolor; and output the fragment color value as required for anyfragment shader program.
 13. A method of claim 9, wherein the fragmentcolor value is output as required for any fragment shader program. 14.The method of claim 1, wherein the channel intensity settings areexecuted by a central processing unit.
 15. The method of claim 1,wherein the lighting modelling method includes a daylight harvestingmodelling system.
 16. The method in claim 1, wherein the lightingmodelling system for displaying real-time dynamic lighting simulationsfor a building.
 17. The method in claim 1, wherein the lightingmodelling system for displaying real-time dynamic lighting simulationsfor horticultural lighting purposes.
 18. The method in claim 16, whereinthe lighting modelling system for displaying real-time dynamic lightingsimulations for a space in a building.
 19. The method in claim 16,wherein the lighting modelling system for displaying real-time dynamiclighting simulations for a greenhouse.
 20. The method in claim 16,wherein the lighting modelling system for displaying real-time dynamiclighting simulations for theatrical lighting.